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SUMMARY 


One of the current goals of research in hypersonic, airbreathing propulsion is access to 
higher Mach numbers. A strong driver of this goal is the desire to integrate a scramjet 
engine into a transatmospheric vehicle airframe in order to improve performance to low 
Earth orbit (LEO) or the performance of a semi-global transport. Airbreathing engines 
such as the scramjet have enhanced performance over non-airbreathing rockets due to the 
removed requirement of carrying oxidizer on board the vehicle. This advantage is reflected 
in enhanced efficiency in the form of specific impulse which translates to lower vehicle gross 
lift off weight when carrying the same size payload to orbit. Additional benefits potentially 
derived from a smaller, lighter vehicle include reduced cost, maintenance, and possibly 
improved reliability. 

There are technical challenges for scramjets flying in the hypervelocity (Mach 10+) 
speed regime. The first issue is residence time of the fuel-air mixture in the combustor as 
the residence time decreases as the freestream velocity of the vehicle increases. The second 
issue is the internal heating loads experienced by the vehicle as velocity increases. In order 
to prevent vehicle materials from exceeding their design temperature limits, active cooling 
is often employed in scramjet designs. This active cooling may have a negative impact to 
efficiency if excess fueling is required. 

The previously discussed problems encountered by scramjets at hypervelocity speeds 
have the potential to be alleviated by engines using premixed, shock-induced combustion 
(PMSIC). In this concept, the fuel is injected on the forebody of the vehicle upstream 
of the inlet. This allows the fuel and air to mix along the forebody instead of in the 
combustor. Once the mixture reaches the inlet, the flow is processed through a shock 
wave of sufficient strength to induce combustion without the use of other ignition devices. 
Because the combustor only needs to support combustion and not mixing its length is 
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reduced considerably, positively impacting the issues of combustor weight and heat loads. 
The burned gases then expand in the vehicle nozzle to produce thrust. 

The goal of the research presented in this dissertation is to understand the flow mecha- 
nisms driving the performance and feasibility of a hypervelocity, shock-induced combustion 
ramjet (shcramjet). This dissertation presents a series of complementary studies simulating 
the physics of a shcramjet vehicle traveling at hypervelocity speeds. Regions of interest 
on the vehicle include the forebody, where injection and mixing must take place with low 
losses and without autoignition occurring. Another region is the inlet, where the presence 
of combustion instabilities must be evaluated. Optimization to maximize performance is 
not within the scope of this research. The design of a specific flight or ground test article 
is also not an objective of this research. 

New knowledge for forebody mixing gained during the course of this research includes 
the comparative performance of the three primary classes of injection at forebody conditions 
relevant to hypervelocity flight. These classes include transverse, ramp, and strut injection. 
A technically defensible means of analyzing the performance of each concept was devised 
and consisted of the direct comparison of the loss-tracking metric as a function of mixing 
efficiency. This provided a valid means to identify injectors that perform the best given any 
particular level of mixing achieved. In doing so, a novel, pathfinding strut injector concept 
was identified as a superior performer relative to the other concepts studied. This concept 
will form the starting point for future studies in forebody injection. 

The presence of autoignition in the mixing regions of the best performing injectors from 
the mixing study were also evaluated. During the course of the study it was found that 
transverse injection, an otherwise high-loss method of injection and mixing, was character- 
ized by extensive regions of autoignition due to shock wave generated by the fuel plume. 
The strut injector had two regions of ignition: a small midstream ignition zone due to the 
presence of shock waves and one at the wall due to the hot boundary layer. Attempts were 
made to mitigate the latter region using both film and transpiration cooling at the wall via 
a parametric study of coolant species and injector conditions. It was found that film cooling 
using either hydrogen or helium provided a significant amount of autoignition mitigation at 
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the wall. 


Shock-induced combustion in a constant area duct was simulated at conditions corre- 
sponding to freestream Mach numbers ranging from 10 to 20 at a constant dynamic pressure 
of 71,820 Pa (1500 psf). It was found that no instabilities were present in the flow in the 
vicinity of the shock wave or elsewhere. At each Mach number simulated, the optimum duct 
length that provided maximum stream thrust potential was extracted. Another study was 
conducted with shock-induced combustion of hydrogen/air ahead of a blunted projectile in a 
quasi-ID framework. As part of this study, the axisymmetric Navier-Stokes equations were 
reformulated in a form valid in the stagnation region of the projectile in order to significantly 
decrease the solution time. A new tool was developed and validated against experimental 
and numerical studies in the literature. It was then applied to parameter studies not yet 
undertaken in the literature due to the feasibility constraints associated with simulating a 
large number of axisymmetric simulations. In doing so, regions of instability were mapped 
out in a large velocity/pressure parameter space. 
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CHAPTER I 


INTRODUCTION 

1.1 Motivation for Current Research 

(What is a Shcramjet and What Are its Challenges? ) 

One of the current goals of research in hypersonic, airbreathing propulsion is access to higher 
Mach numbers. A strong driver for attaining this goal is the desire to integrate a scramjet 
engine 1 into a transatmospheric vehicle airframe in order to improve single stage or two 
stage performance to low Earth orbit [1], Airbreathing engines such as the scramjet enjoy 
enhanced performance over non-airbreathing rockets due to the eliminated requirement of 
carrying oxidizer on board the vehicle. This advantage is reflected in enhanced efficiency 
(i.e. I sp ) which translates to lower vehicle gross lift off weight when carrying the same sized 
payload to orbit [2]. Additional benefits potentially derived from a smaller, lighter vehicle 
include reduced cost, maintenance, and possibly improved reusability. 

Scramjets are operationally limited to the atmospheric stage of flight due to the re- 
quirement that the atmosphere provide the oxygen necessary for combustion. Transition to 
rocket mode to complete the flight to orbit occurs at the point in the trajectory where the 
air is too thin for scramjet operation or earlier in the trajectory due to other constraints. 
Because the subsequent pull-up stage is of lower efficiency, the staging Mach number is 
critical in defining the overall efficiency of the vehicle; the higher the staging Mach number 
the more efficient the overall vehicle. Figure 1 demonstrates the difference in staging flight 
velocity for a vehicle at Mach 10 versus one at Mach 15. The higher first stage maximum 
Mach number achieves 38% of the required kinetic energy to LEO versus 15% for the lower 
first stage maximum Mach number. 

1 Scramjet engines are characterized by a flowpath with ideally no moving parts and air compression 
accomplished by shock waves produced by the forebody and inlet of the vehicle. This will be discussed in 
more detail later. 
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Figure 1: Staging altitudes and velocities for a vehicle with a first stage maximum speed 
of Mach 10 versus a vehicle with a maximum first stage speed of Mach 15. 

There are technical challenges that scramjets face when attempting to fly in the hyper- 
velocity 2 speed regime. The first issue is related to the residence time of the fuel-air mixture 
in the combustor. The residence time of the fuel within the combustor (r r « L c /V 0 c ) rapidly 
decreases as the freestream velocity of the vehicle increases. This severely limits the time 
available for fuel injection, fuel-air mixing, and burning within the combustor. While the 
combustor could be lengthened to compensate for decreasing residence time, doing so will 
create additional weight burdens on the vehicle as the combustor tends to be a significant 
weight component of the engine. 

The second issue is related to the heating loads experienced by the vehicle. As velocity 
increases, the heating loads experienced by the external and internal surfaces of the vehicle 
become more harsh. In order to prevent vehicle materials from exceeding their design 
temperature limits, active cooling is often employed in scramjet designs. Primary methods 
of active cooling are 

1. Passing cryogenic hydrogen though channels below the flowpath surface 

2. Injecting a film of hydrogen gas at the surface in the internal flowpath (so-called “film 
cooling” ) . 

2 The term hypervelocity typically refers to Mach numbers above 10, which corresponds to velocities in 
excess of 3 km/s on 1000 to 2000 psf trajectories. In this velocity regime, stagnated regions will begin to 
exhibit oxygen dissociation [3]. 
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Both of these cooling methods provide lower surface temperatures with the first method 
preheating the fuel before it enters the combustor. As the vehicle velocity increases past a 
critical limit, however, a fuel system supplying near-stoichiometric hydrogen to the flowpath 
is no longer sufficient for cooling and an increase in fuel flow rate is necessary to keep the 
surfaces cool [2]. This causes two negative impacts. First, the excess fuel is wasted for lack 
of oxidizer with which to react. Second, vehicle specific impulse is negatively impacted by 
the increased fuel mass flow rate, reducing engine efficiency. 

Attention is now given to the current concepts occupying the hypersonic and hyperve- 
locity regimes of flight. The first is the supersonic combustion ramjet (i.e. scramjet), shown 
notionally in Figure 2. Because of the highly integrated nature of hypersonic propulsion, 
the term scramjet is typically used to describe not only the engine but the vehicle itself. 
However, here the discussion is limited to the engine only. The typical scramjet concept is a 
single, continuous flowpath ideally with no moving parts (e.g. no mechanical compressors or 
turbines) . The incoming airflow is compressed on the forebody of the vehicle via a system of 
oblique shock waves that, for an ideal case, terminate at the shoulder of the inlet 3 . Once the 
compressed air passes into the fully enclosed flowpath, some concepts require the use of an 
isolator section that is characterized by a shock-train interacting with a growing boundary 
layer (e.g. a “dual-mode” scramjet). The isolator section provides the necessary time for 
the pressure rise demanded by the combustor so that inlet unstart does not occur. Once the 
air enters the combustor fuel is injected where it is then mixed with the air stream. Com- 
bustion occurs by igniting the mixture, for example with a silane- hydrogen pilot stream or 
through autoignition if the ambient conditions allow. The burned gas mixture then expands 
through the nozzle to create thrust. 

The previously discussed problems encountered by scramjets at hypervelocity speeds 
may be alleviated by engines using premixed, shock-induced combustion (PMSIC), shown 
notionally in Figure 3. In this concept, the fuel is injected on the forebody of the vehicle 
upstream of the inlet. This allows the fuel and air to mix along the forebody instead of in 

3 The shoulder is the point on the vehicle forebody where the flow experiences a local expansion. Two 
dimensional designs typically aim to terminate the cowl shock wave on the shoulder to minimize shock 
structures in the internal flowpath. 
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Isolator shock train 
Figure 2: A notional scramjet 


the combustor. Once the mixture reaches the inlet, the flow is processed through a shock 
wave of sufficient strength to induce combustion without the use of other ignition devices. 
Because the combustor only needs to support combustion and not injection and mixing, its 
length is reduced considerably, positively impacting the issues of combustor weight and heat 
loads. The burned gases then expand on the vehicle nozzle to produce thrust. Engines using 
PMSIC are called either shock-induced combustion ramjets (e.g. shcramjets) or oblique 
detonation wave engines (e.g. ODWEs, or wave engines). The distinction between the two 
is in the distance between the combustion front and the shock front inducing it: shcramjets 
typically have a measurable distance between the two while wave engines do not (i.e. the 
inducing shock wave is a detonation wave). The remainder of this dissertation will use the 
term shcramjet to refer to an engine using PMSIC. 


Fuel Injectors 



Cowl Lip Combustion Zone 
Figure 3: A notional shcramjet 

Although it addresses many of the problems faced by scramjets, the shcramjet concept 
encounters technical challenges of its own. These challenges are illustrated in Figure 4 and 
summarized as follows: 

Nonuniform forebody mixing Fuel injected on the forebody of the vehicle must be 
able to mix effectively into a cross stream of very high Mach number (e.g. Mach 7-8 for 
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a Mach 12 vehicle). The effect of high Mach numbers on vorticity and mixing of the fuel 
jet is not well understood. Additionally, uniform mixing is desired in order to allow even 
burning in the combustor of the shcramjet. 

Boundary layer autoignition As the injected fuel enters the hot forebody boundary 
layer there is a risk of autoignition of the resulting fuel-air mixture near the body. This 
would result in the generation of severe pressure drag, high wall heating, and significant 
fuel losses. Therefore, the prevention or mitigation of autoignition on the forebody is of 
high importance. 

Shock wave instability It is possible for combustion instabilities to develop in the reac- 
tion zone due to interaction between the shock wave and the induced reaction front. This 
has been observed in fundamental studies of both blunt body and wedge-supported shock- 
induced combustion. Depending on the conditions of the oncoming flow, the system may 
either be completely stable or unstable with possibly large-amplitude behaviors. This can 
occur in a perfectly mixed system although non-uniformities in mixing may contribute to 
instability as well. Consequences of these instabilities may range from performance loss to 
engine unstart. 

Nonuniform mixing and/or spillage of fuei/air mixture 



Figure 4: Technical challenges for a notional shcramjet. 


The goal of this research is to evaluate and understand the mechanisms driving the 
performance of a hypervelocity shcramjet vehicle. This research presents a series of comple- 
mentary studies that simulate the physics of a notional shcramjet vehicle traveling at and 
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above Mach 12. Regions of interest on the vehicle include the forebody, where injection and 
mixing must take place without autoignition, and the inlet-combuster, where combustion 
instabilities must be minimized to ensure predictable and uniform performance. These re- 
gions of interest on the vehicle will be studied individually with inflow boundary conditions 
traceable to hypervelocity freestream flight conditions. The objective of this study is a 
physics-based understanding of the performance drivers of the shcramjet system; however, 
black-box optimization to maximize performance is not within the scope of this research. At 
the conclusion of this research, understanding of the physics driving hypervelocity vehicle 
performance will be enhanced, thus allowing an analyst to effectively design a shcramjet 
flowpath. 

1.2 Background of PMSIC Technology 

This section provides orientation to the high-level background of PMSIC as related to 
the identification of cycle concepts using the technology and their benefits over the state- 
of-the-art of the time as well as today. While high-level system studies raise important 
questions (e.g. How can fuel be injected while minimizing losses? How may fuel be prevented 
from burning in the vicinity of the injector?), it will be left to subsequent chapters of this 
dissertation to provide background and study to them. 

Initial cycle analyses of the detonative ramjet originated in the 1950’s and were limited to 
top-level parameter studies assuming ideal mixing and burning. Dunlap, et al. [4] evaluated 
the benefits of an airbreathing engine using wedge-stabilized, detonative combustion versus 
other, more conventional cycles. Their concept of such an engine is shown in Figure 5. Their 
analysis assumed perfectly mixed fuel and air upstream of the detonation location with no 
losses present except for those caused by the detonation wave. The authors brought up 
the possibility of an unstable detonation wave, concluding that a wedge centerbody would 
in all likelihood provide the necessary stabilization. They also bring up the possibility for 
fuel burning in the boundary layer in the vicinity of the fuel injectors. This issue persists 
today and motivates the current research. It was realized in their analysis, assuming ideal 
mixing and a fixed ignition temperature, that a hydrogen-fueled concept can realize similar 
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performance to a conventional ramjet at higher Mach numbers. While the study suffered 
from drawbacks such as not understanding the autoignition properties of the fuels used, it 
was an important initial data point in the literature. 



Figure 5: A reproduced schematic of the concept analyzed by Dunlap, et al. [4] is shown. 

Soon after the Dunlap study, Sargent and Gross [5] published a study analyzing the 
detonation wave ramjet in further detail using hydrocarbon fuel. The study was derived 
from experimental results in plane and oblique detonation, reported by Gross [6, 7]. As 
with the Dunlap cycle, the engine concept was assumed to compress the incoming air 
stream to a target supersonic velocity subject to varying inlet kinetic energy efficiencies. A 
representation of the engine analyzed by Sargent and Gross is shown in Figure 6. Assuming 
no losses from injection, uniform mixing, and complete combustion, the I sp , SFC, thermal 
efficiency, and exit to inlet area ratio were determined as a function of flight and detonation 
Mach numbers, inlet efficiency, and nozzle efficiency. The results indicated that compared 
to the subsonic burning ramjet, the detonation wave engine has lower fuel consumption 
at the cost of thermodynamic performance (which is a minor difference for some cases). 
The report also identified other salient points relating to the ODWE, for example shorter 
combustor lengths (due to reduced mixing and reaction length), lighter structure (due to 
lower static pressure), and lower heat loads (due to shorter lengths). The problem of fuel 
injection was also identified, an issue that persists today and motivates the current research. 
The need for variable geometry was also identified as a critical need. 

Rubins and Rhodes [8] conducted experiments in premixed, shock-induced combustion 
at the Arnold Engineering Development Center in 1963 that were important pathfinding 
studies in the feasibility of the technology. Rather than relying on the creation of a detona- 
tion wave for combustion, an oblique shock was used to raise the temperature of a premixed 
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Flow 



Figure 6: A reproduced schematic of the concept analyzed by Sargent and Gross [5] is 
shown. 

fuel-air stream above the mixture’s autoignition temperature. In this sense, what they were 
proposing was not a “wave” engine, but a more general shock-induced combustion process. 
A flowpath using this technology is shown in Figure 7. The process shown differs from an 
oblique detonation in the sense that the reaction zone is no longer tightly coupled with the 
shock wave, the only purpose of which is to increase the temperature of the fuel-air mixture. 



Figure 7: A reproduced schematic of the concept proposed by Rubins and Rhodes [8] is 
shown. 

In order to test the concept, Rubins and Rhodes created an experiment that mixed 
fuel and air before the mixture was processed through an oblique shock supported by a 
28° wedge. A schematic of the experiment is shown in Figure 8 with the upstream flow 
provided by a combustion-heated Mach 3 tunnel. The experiments showed that it is possible 
to control chemical reaction through the use of shock waves and demonstrated predictable 
behavior in the combustion process compared to theory. Similar to previous studies, Rubins 
and Rhodes identified the issue of fuel injection as a critical one moving forward. 

Rubins and Bauer [9] made cycle calculations to evaluate the benefits of preinjection as 
applied to scramjet flowpaths-as sketched in Figure 9-in the Mach 9 to 22 range. Their 
study was distinct in that the concept studied had fuel injected on the forebody of the 
vehicle during the compression cycle. While one proposal for ignition of the premixed 
fuel-air stream was via shock waves, it was only one of a handful of other options, for 
example the more conventional pilot flame/jet. The study was completely analytical in 
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Figure 8: A reproduced schematic of the PMSIC experiment conducted by by Rubins and 
Rhodes [8] is shown. 

nature and identified key benefits of forebody injection, including improved fuel specific 
impulse due to lowered forebody wall friction (for a lighter-than-air fuel) and improved 
pressure recovery. Methods to limit forebody fuel ignition were also discussed, for example 
parallel fuel injection and the use of him cooling. 

Shock 



Figure 9: A reproduced schematic of the concept proposed by Rubins and Bauer [9] is 
shown. 

Through the 1970’s, PMSIC research generally concentrated on unit studies of shock- 
induced combustion ahead of blunted projectiles and wedges [10, 11, 12] without much 
direct maturation of engine concepts using PMSIC. These intervening studies did provide 
important steps in understanding the physics of the phenomenon, as will be discussed in 
later chapters. As for the lack of further concepts using PMSIC, research focus most likely 
drifted away from the technology due to near-term challenges (e.g. forebody auto-ignition) 
and the fact that scramjet flight at low- to mid- hypersonic (i.e. 5 - 10) Mach numbers was 
and is realizable without premixing. 

The period from the late 1970’s through the early 1990’s brought important studies 
toward the development of PMSIC technology. First, there were a number of systems 
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studies [13, 14, 15, 2] that brought the concept of a transatmospheric vehicle using PMSIC 
engines into greater focus and showed their potential benefits over “traditional” scramjet 
engines. These detailed studies demonstrated an important pull case for the technology 
and remain important references for the quantitative systems benefit of PMSIC. Second, 
the advent of the National AeroSpace Plane (NASP) concept was a milestone in the national 
drive toward developing an orbital access vehicle using airbreathing engines. The orbital 
access use case is an important one for PMSIC technology. Third, through the same time 
period, computational fluid dynamics experienced rapid improvements in capability and 
some of the first numerical simulations on PMSIC and its application to an engine were 
conducted [16, 17, 18, 19]. Many of these initial studies marked the beginning of long-lasting 
and geographically diverse research programs into the technology. Finally, experiments in 
PMSIC found some resurgence [20, 21] and further demonstrated the feasibility of the 
technology. 

With the systems benefits of shcramjet engines having been established, research in the 
physics of PMSIC has experienced renewed focus in the new millennium. Fundamental 
studies in blunt-body PMSIC have maintained attention in Japan and wedge-supported 
PMSIC has continued to be studied in the United States and elsewhere. In addition, 
research programs addressing systems issues related to PMSIC, for example fuel injection, 
has received attention by researchers in the United States, Canada, and Europe. While the 
background of this research will be outlined in future chapters, the fact that strong research 
programs exist to develop the technology for use in hypervelocity engines lends optimism 
for the future of the shcramjet concept. 

1.3 Knowledge Gaps in the Literature 

While fundamental advances have been made in advancing the state of the art in PMSIC and 
its application to shcramjet engines, the literature currently has several gaps that are open 
avenues for pathfinding research. Selected knowledge gaps are categorically summarized as 
follows: 


10 



Goal: Efficient (i.e. high gain, low loss) injection and mixing performance with 

stable combustion 


Figure 10: Studies involved in the present work and knowledge to be gained. 

Forebody Injection/Mixing: Few concepts have been studied in depth for fore- 
body fuel injection in the hypervelocity regime. Recent studies have concerned the 
forebody injection of fuel using either a slot [22] or modified ramps [23] . While a num- 
ber of studies from the same research group have tackled the latter geometry in terms 
of characterizing its performance, any geometrical modification has been through per- 
turbation of the same design. No studies exist in the literature comparing several 
injection concepts at the same conditions in order to make a direct comparison. 

There is also a lack of proper characterization of injector performance in the literature. 
In the general literature of fuel injection and mixing, many studies neglect to make any 
characterization of the losses incurred by the injector concept. The typical means of 
measuring the efficacy of a particular concept is to limit the analysis to positive metrics 
such as mixing efficiency, or to a more coarse extent, fuel penetration or injected 
plume area. By ignoring losses, it is not possible to truly gain an understanding 
of the relative performance of different injectors, as will be discussed in Chapter 2. 
Of the studies that do attempt an analysis of losses, many use the common metric 
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of total pressure loss. While this metric is appropriate for determining the losses 
incurred during a general compression process, it is not appropriate where the mixing 
of gases is occurring (this will also be discussed in Chapter 2). Finally, of the studies 
that do make a proper characterization of losses, there is little clarity as to how to 
properly compare positive metrics to loss metrics. This is because mixing itself is a 
loss- producing process. 

Autoignition Mitigation: The tendency for fuel injected on a vehicle forebody to 
ignite when exposed to the hot boundary layer region has long been identified as a po- 
tential problem for such concepts, but few options for mitigation exist in the literature 
as pertains to hypervelocity propulsion. To the author’s knowledge, the only concept 
for autoignition mitigation for a hypervelocity forebody is the slot injection of nitrogen 
at the base of a cantilevered ramp injector [24]. While it demonstrates promising re- 
sults, other mitigation strategies may provide similar or better performance for lower 
amounts of injected transpirant. 

Shock Instability: Shock instability at the combustor entrance has two potential 
sources. The first is due to the interactions and dynamics between the shock wave 
and the resulting reacting front. This phenomenon has been studied to a degree 
using blunt-body problems, but there is no comprehensive parametric study on what 
the effect of composition, temperature, pressure, and velocity has on the resulting 
dynamics. This is due to the logistics of the problem, where most researchers create 
axisymmetric computational grids with convergence requirements much higher than 
that which is needed for a rapid parametric study. The second instability is due to 
an unstable solution to the oblique shock problem. This issue has been looked at to 
some degree by researchers, but the effect of this instability of fuel composition and 
wedge bluntness is not well characterized (although Matsuo conducted a small study 
of PMSIC on a blunted wedge for two bluntness radii). 
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1.4 Objectives and Contributions 


It is the goal of this dissertation to address the knowledge gaps outlined in the previous 
section. The same categories of knowledge gaps will compose the major sections of this 
thesis. The first objective of this dissertation will be to compare injection methodologies at 
the same flow conditions to determine relative benefits of each. Note that the objective of 
this study is not to provide optimized solutions to the problem of forebody injection and 
mixing, but rather to conduct comparative studies to identify methods that outperform 
others. In addition, quantitative benefit and loss metrics will be compared directly in order 
to determine the value of one concept over another. 

The second objective will be to address the issue of autoignition in the vicinity of the fuel 
injectors. In order to reduce the number of injectors considered, only the best performing 
injectors from the comparative study will be evaluated. Film and transpiration cooling 
will be evaluated for their ability to mitigate autoignition and a parametric study will be 
conducted to understand the effect of different coolant flow variables. 

Finally, the issues of shock/reaction front stability will be investigated. The axisymmet- 
ric Navier-Stokes equations will be reformulated in a form that lends itself to rapid analysis 
with minimal loss of fidelity. The shock/reacting front dynamics will be characterized with 
respect to the important parameters affecting the flow. Furthermore, it will be determined 
if instabilities occur in a constant area duct representative of a notional combustor. 

The specific, unique contributions of this research that advance the body of work in this 
flight regime include: 

1. Transverse, ramp, and strut injection concepts simulated at the same conditions to 
determine... 

(a) ...mixing efficiency and thrust potential as a function of distance from injector 
leading edge. 

(b) ...the relative performance of each concept using the ratio of mixing efficiency 
and thrust potential as a function of distance from the injector leading edge. 
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(c) ...where autoignition is present in the vicinity of best-performing injectors and 
the mechanisms by which autoignition occurs. 

2. How best to mitigate autoignition in the vicinity of the wall using film cooling and 
wall transpiration, including... 

(a) ...the best selection of transpiration locations relative to the injector. 

(b) ...the proper selection and effect of cooling parameters such as wall porosity, mass 
flow rate, species, and coolant temperature. 

3. Determination of the optimum notional combustor length for different points on a 
notional airbreathing trajectory as well as evaluating if there are instabilities present 
in the internal flow. 

4. Improved predictions of the regimes of shock/reacting front stability by enabling a 
rapid simulation framework by the re- formulation of the Navier Stokes equations in a 
quasi- ID form. 

These objectives will be accompanied by numerical studies to validate the specific physics 
controlling the behaviors of each type of flow (e.g. mixing, combustion, unsteady dynamics) 
as simulated by the CFD code chosen for the analysis. 

1.5 Dissertation Preview 

Chapter 2 concerns the injection and mixing of fuel with air at conditions relevant to 
a hypervelocity vehicle on a 71,820 Pa (1500 psf) trajectory at the Mach 12 trajectory 
point. Injection is studied in nonreacting simulations among the major categories of 
injection: flush-wall, ramp, and strut injection. 

Chapter 3 concerns the mitigation of autoignition on a notional shcramjet vehi- 
cle forebody. For this study reacting simulations of film and transpiration cooling 
concepts are evaluated as a means to cool the boundary layer. 

Chapter 4 concerns the stability of the shock and reaction front located at the 
entrance of a notional shcramjet combustor entrance. Optimum combustor lengths 
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will also be obtained as a function of point on the trajectory using stream thrust 
potential. 

Chapter 5 concerns the development of a tool to rapidly analyze the stability char- 
acteristics of a premixed, shock-induced combustion shock and reaction front ahead 
of a supersonic, blunted projectile. 

Chapter 6 provides a summary and conclusion to the research presented herein, as 
well as avenues for future research. 
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CHAPTER II 


FOREBODY INJECTION AND MIXING 


Much of the literature on fuel injection research for scramjet combustion has been at Mach 
numbers representative of mid-range hypersonic flight (Mach 5-8). Of these, most numer- 
ical and experimental studies have been conducted at conditions relevant to the combustor 
as typical fuel injection schemes begin at the transition from the isolator 1 to the com- 
bustor. Therefore, most studies have taken place around Mach 2 - 3 in duct flow. In 
this speed regime, a significant number of experiments have evaluated perpendicular wall 
injection, [25], near-parallel injection concepts (e.g. swept ramp injectors) [26, 27], transpi- 
ration/effusive injection cooling for transverse injection [28], and the effect of fuel injector 
nozzle geometry on mixing characteristics [29, 30, 31]. The use of cavities for injection and 
flameholding has also been a subject of past research [32], The use of oblique shock waves 
to enhance mixing has been proposed, tested, and used in practice [33, 34], 

Other studies have investigated the use of pre- injection on the vehicle forebody in or- 
der to allow more time for mixing with the goal of reducing combustor length. Rubins 
and Bauer [9] described an early analytic treatment of pre-injection with application to 
scramjet engines from Mach 9 to 22 flight conditions and concluded that premixing with 
hydrogen has the ability to improve inlet pressure recovery, reduce friction, and shorten the 
vehicle combustor. A numerical study of slot-injection into a Mach 6.3 flow was conducted 
by Gonzalez [22], investigating the parametric effect of fuel temperature and pressure on 
the potential for ignition in the boundary layer. Recent advances in pre-injection for the 
shcramjet forebody with hydrogen and hydrocarbon fuel considered the use of cantilevered 
injectors with slot cooling to prevent boundary layer autoignition [35, 24, 36, 23, 37, 38]. 

A schematic of the flowfield about a transverse, underexpended fuel injector is shown in 

1 The isolator section of a scramjet is a component to protect the inlet of the vehicle from the pressure 
rise in the combustor for low supersonic Mach numbers. Without the isolator, engine unstart may occur at 
some conditions. 
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Figure 11: A schematic of the flow in the vicinity of an underexpanded transverse fuel 
injector in supersonic flow is shown. 

Figure 11. The penetration of a transverse fuel jet into a crossflow is largely governed by 
the dynamic pressure ratio, J, of the fuel jet to the crossflow (e.g. relative force of blowing) 
[25], the geometry of the jet itself (e.g. disturbance to the crossflow) [30, 39], and the angle 
of injection. The dynamic pressure ratio is bounded by maximum fuel system pressure 
and the desire to produce an underexpanded plume. An underexpanded plume is desirable 
because it prevents losses due to shock structures forming within the jet. The geometry of 
the jet is governed by the injector nozzle and its exit shape with more complex geometries 
representing a more difficult problem in fabrication and survivability. The angle of injection 
represents a design tradeoff between high-penetration/high- loss normal injection and low- 
penetration axial injection that imparts a beneficial thrust to the vehicle due to the fuel jet 
itself. 

Identifying an injector that best penetrates and promotes mixing of fuel into the main 
flow is insufficient when considering how the injected plume is to subsequently mix with 
its surroundings. Molecular mixing is governed by diffusion, and how quickly fuel and air 
diffuse through a unit area of interface is largely out of the designer’s control. Therefore, 
increasing the mixing rate through the fuel-air interface is accomplished by stretching the 
interface using vorticity in the flow. Vortical stirring structures can be induced via geometry 


17 


placed in the flowfield as is the case with the swept ramp concept [26] and by flow interaction 
phenomena such as that between a fuel jet passing through a shockwave [34], Care should 
be taken in designing such an interacting flowfield, however, as stirring phenomena can 
work against the designer. For example, some character of shock waves can turn fuel back 
toward the wall [35]. 

Once a fuel injector is identified that has excellent penetration and mixing characteris- 
tics, the designer is once again presented with a challenge in the form of losses produced 
by the injector concept. For example, there may exist an injector concept that is able to 
achieve perfect mixing but at the cost of depleting the main flow of the momentum upon 
which the vehicle depends for producing thrust. There are a number of loss sources that 
must be taken into account when designing and implementing a forebody injector. For 
example, separation zones are likely to develop that may entrain the fuel from the injector 
and enhance local combustion. In particular this is observed in flush- wall injection [40]. 
This is a potential source of drag from combustion on the forebody as well as fuel loss. 
Combustion in the vicinity of the injector also has the potential for increasing local heating 
rates. Another source of drag is autoignition due to the fuel entering the hot boundary 
layer, igniting, and pressurizing the forebody. This is a significant source of drag, fuel loss, 
and intense heating at the wall. This issue can be mitigated with film or transpiration 
cooling; however, care should be taken in this approach due to the negative impact on the 
specific impulse with excess fueling. Other sources of drag include intrusive geometries and 
unpressurized, rearward facing steps. Finally, injection schemes that deviate from axially 
directed injection suffer increasing momentum losses as the angle of injection approaches 90 
degrees. Some of these factors may actually be a benefit in the context of internal injection 
(e.g. rapid autoignition) but are detrimental for forebody injection. 

The purpose of this chapter is to report on the study of forebody fuel-air mixing on 
a notional forebody at conditions representative of hypervelocity flight. Hydrogen fuel 
injection was simulated on flat plate with the crossflow conditions representing conditions 
of a Mach 12 point on a 71,820 Pa (1500 psf) trajectory after being shocked through two 
forebody shocks. Flush wall, ramp, and strut injection were simulated numerically and 
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compared using a novel methodology for reducing the mixing and loss metrics extracted 
from the flow. 

2.1 Numerical Framework and Approach 

While a great deal of effort could have been spent to develop a new CFD solver for the 
work contained in this dissertation, it was decided instead to leverage the capabilities of an 
existing numerical code. This is because the main objective of the mixing and autoignition 
studies is physical investigation, not tool development. Therefore, the numerical studies 
in Chapters 2 and 3 were carried out using the Viscous Upwind aLgorithm for Complex 
flow ANalysis (VULCAN) [41] CFD code, developed and maintained at NASA Langley 
Research Center. VULCAN uses a finite- volume, cell-centered scheme for solving calorically- 
or thermally-perfect flows on a structured grid. While some CFD codes permit the use 
of unstructured grids, beneficial when creating topologies for complicated geometries, the 
structured grid capabilities of VULCAN allow for generally higher accuracy for fewer grid 
points. This is especially true when evaluating viscous regions which are important in this 
dissertation. 

VULCAN began development in 1992 under the name Langley Algorithm for Research 
in Chemical Kinetics (LARCK) as a replacement of the SPARK code [42] , motivated by the 
desire to study complex geometries and incorporating upwind schemes into the simulation 
capability 2 . Instead of modifying SPARK to suit evolving needs for high speed propulsion 
applications, LARCK began development without a heritage codebase. Because turbulence- 
chemistry interaction problems were of interest at the time, LARCK also included models 
for such phenomena. Early in its development, LARCK was successfully validated against 
multidimensional flat plate and ramp flows as well as used in the simulation of integrated 
scramjet flowpaths [43]. In 1996, LARCK was upgraded with new space marching capa- 
bilities and turbulence wall functions under contract from the Air Force and it was at this 
time the code was renamed to VULCAN 2 . In the late 1990s, an effort was undertaken 
to parallelize the code. Some examples of validation studies undertaken with VULCAN 

2 Jeffrey A. White, personal communication 
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include Mach 3 flow over a cavity flameholder [44], Mach 1.8 coaxial jet flow [45], Mach 2 
vitiated flow in a facility nozzle [46], a Mach 4 sidewall compression inlet [47], and a Mach 
6 streamline traced inlet [48]. 

The VULCAN code is designed for and excels at the computation of high speed, exother- 
mic flows characteristic of scramjet engines. The algorithm is valid from roughly Mach 0.1 
to any speed regime that would cause a thermochemical breakdown in the flow. For exam- 
ple, while at the time of this writing vibrational nonequilibrium is being integrated into the 
code, it is currently unsuitable for the simulation of very high speed endothermic flows such 
as those studied in entry, descent, and landing applications. Because of its development 
and specialization in the simulation of high speed, exothermic flows, VULCAN has entered 
wide use for both fundamental hypersonic research and in successful test programs such as 
Hyper-X (X-43A/B/C), X-51, and HIFiRE. 

The flow for this study was assumed to be viscous and thermally perfect. Non-reacting 
simulations were carried out in the mixing studies of this chapter and reacting simulations 
and the effect of autoignition were addressed in Chapter 3. Selected capabilities of the 
VULCAN CFD code are summarized in Table 1. In this study the flow was integrated 
spatially using the LDFSS scheme [49] to improve simulation robustness and discontinuity 
resolution in the flow field. The MUSCL interpolation scheme was chosen to be third order 
upwind-biased to reduce truncation error and the use of a van Leer flux limiter assured 
that the solution is total variation diminishing (TVD). Unless otherwise noted, the flow 
was integrated temporally using the DAF scheme [50] with local time stepping and an 
increasing Courant-Friedrichs-Lewy (CFL) number schedule. The use of DAF allows for 
more stable behavior during solution although ILU increases the CFL number limit at the 
cost of significantly increased storage requirements. 

Unless otherwise stated, the simulated flow was assumed to be fully turbulent. VUL- 
CAN is equipped for Reynolds- averaged simulation (RAS) using the one equation turbulence 
model of Spalart or two equation models of Menter (with and without shear stress exten- 
sion) or Wilcox. Recent VULCAN developments include simulating turbulence using large 
eddy simulation (LES) and hybrid RAS/LES. While the LES and hybrid methods both 
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Table 1: A capability summary of the VULCAN CFD code 
Spatial integration schemes 
Parabolic (space marching) 

Elliptic 

Inviscid flux reconstruction schemes 

None 

Roe flux differencing 
Local Lax Friedrichs 
van Leer flux vector splitting 
Harten, Lax, van Leer, and Contact flux splitting 
Edwards low dissipation flux splitting 
Edwards low Mach number preconditioned flux splitting 
Flux limiters 
None 

Minmod (TVD) 
van Leer (TVD) 
van Albada 
Smooth 
Koren 

Minmod (ENO) 
van Leer (ENO) 

Turbulence models 

Laminar 

One-equation incompressible Spalart 

Two-equation Menter blended Wilcox 1988 k-w/ Jones Launder k-e 
Two-equation shear-stress transport Menter variant 
Two-equation 1998 Wilcox k -uj 
Two-equation 2006 Wilcox k-w 

Explicit Algebraic Stress model of Rurnsey and Gatski (2003) 
Constant coefficient Smagorinsky sub-grid scale model 
with van Driest wall damping 

Chemistry modeling 

Frozen 

Finite Rate (Arrhenius with or without Troe falloff parameters) 
Computer Automated Reduced Mechanism (CARM) 

CARM with In-Situ Automated Tabulation (IS AT) 
Point-implicit chemistry treatment 
Analytical Jacobian 
Numerical Jacobian 

Temporal integration schemes (local or global) 

Diagonal Approximate Factorialization (DAF) 

Implicit LU-Factorialization (ILU) 

Other options 

Multi-grid cycling 

2D, Axisymmetric, or 3D domains 
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provide improved accuracy and physical reproduction, for this dissertation it was decided 
to treat turbulence using RAS in order to take advantage of the significantly improved sim- 
ulation time of the two-equation methods. Therefore, turbulence calculations in this study 
use the two-equation Menter k-c j model [51] with a turbulent Prandtl number (governing 
turbulent energy diffusion) of 0.9, turbulent Schmidt number (governing turbulent mass 
diffusion) of 1.0, and a turbulence intensity of 0.1. The turbulence parameters were held 
constant for each simulation to ensure consistency of comparisons made between cases. For 
wall-bounded flows, VULCAN has options for either solving-to-the-wall or employing wall 
matching functions (WMF). It was decided to use the wall matching function of Wilcox 
[52] in order to relax the grid spacing requirements at the wall. Without WMF, the first 
point in the topology typically needs to be located in the viscous sublayer at a y + value 
between 1.0 and 2.0 while the use of matching functions relaxes this requirement to some- 
where between 20.0 and 200.0 [53] with lower values corresponding to better prediction. For 
y + values under 10.0, VULCAN uses a blending scheme to transition from wall matching 
to a solve-to-the-wall method. 

2. 2 Validation 

The process of validating a tool is an important step to ensure that a particular code 
or model is accurately reflecting the physics of the modeled regime. The Department of 
Defense [54] defines validation as 

“Validation: The process of determining the degree to which a model is an 
accurate representation of the real world from the from the perspective of the 
intended uses of the model. ” 

In other words, a tool has been validated for use in a particular regime if simulations carried 
out with the model at the same conditions accurately match reality through the use of 
validation data (e.g. experimental data). Validation should not be confused with verification 
(i.e. establishing that the programmed model matches the conceptual/theoretical model) or 
qualification (i.e. establishing that the conceptual/theoretical model matches reality with 
a defined level of fidelity) . 
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The computational fluid dynamics software of choice in this study was VULCAN, a 
NASA Langley Research Center code specializing in the simulation of supersonic, exother- 
mic flows. More information on the code was given in Section 2.1. The ability of VULCAN 
to properly predict the physics of fuel injection into a quiescent or flowing reservoir was 
verified by comparing simulated results against theory and experiment at similar conditions. 
The validation questions are summarized as: 

Can VULCAN... 

• ...successfully predict the Mach disk height, Mach disk width, and wavelength of the 
jet plume when injected into quiescent air at various jet pressure ratios? 

• ...successfully predict the Mach disk height of a jet plume injected into a supersonic 
crossflow? 

• ...successfully predict the mixing plume characteristics for axial injection into a su- 
personic coflow? 

These questions were chosen because the important controlling physics for the injection 
and mixing studies to follow include the displacement of injected fuel jets from the wall as 
well as the effect of vorticity on the mixing of fuel with air. If the answer to each of these 
validation questions is “yes,” it can be reasonably concluded that VULCAN can predict 
the physics of injection and mixing in the present study assuming similar flow features are 
present and/or the conditions are comparable. 

2.2.1 Transverse Injection into a Quiescent Medium 

The ability of VULCAN to simulate the physics of sonic jet injection into a quiescent 
medium will now be summarized. This problem has been studied in the past by many 
researchers numerically, experimentally, and theoretically and is therefore a good validation 
problem to give confidence to VULCAN’s capability in predicting such physics. There is 
similarity between the current validation problem and the transverse injection and mixing 
study contained in this thesis in that a diatomic gas is being injected at sonic conditions, 
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Figure 12: The grid used for injection validation is shown with every fourth grid point 
drawn. For cases with crossflow, the flow originates from the upper right. 

although the speed of the medium is different between the validation case. The primary 
purpose of this study was to evaluate VULCAN’s ability to predict the structure of the 
barrel shock of the injected jet. 

A cutaway of the topological grid used in this study is shown in Figure 12 with every 
fourth grid point drawn. The grid had 2.5 million grid cells and was clustered near the 
injector, the wall, and the symmetry centerline. A detail view of the circular injector in 
Figure 13 shows the grid topology of the injector, which consisted of a C-grid wrapped 
around a H-grid. Studies involving jet injection into a quiescent medium typically use an 
axisymmetric grid that extends radially from the plane of the jet exit and axially down- 
stream of the jet centerline. It was decided, however, to use the same grid for both cases of 
a sonic jet issuing into a quiescent and supersonic crossflow. This also reduced the presence 
of instabilities typically found in axisymmetric simulations of jets issuing in to quiescent 
mediums. 

Experimental evaluations of sonic jets containing 95% N2/5% NO by volume into a 
quiescent chamber were done by Wilkes et al. [55] at NASA Langley Research Center. The 
Nitrogen jets were seeded with NO in order to provide visualization capability by planar 
laser-induced fluorescence of nitric oxide (NO-PLIF). The visualization method used was 
able to capture clear images of the structure of each underexpanded jet into a quiescent 
chamber. Multiple jet conditions were tested that provided data of jet primary wavelength 
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Figure 13: Injector detail for the validation grid with all grid points drawn. The injector 
topology consists of an H-grid near the center of the injector with a C-grid wrapped around 
it. 
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as well as Mach disk width and standoff distances as a function of the pressure ratio between 
the jet and the chamber. Therefore, these experimental results are a good source of high- 
quality data against which to compare the present numerical results. Figure 14 shows an 
annotated bisection of a jet from a converged VULCAN simulation with a total pressure 
ratio of 20. 

The empirical and theoretical results of Love, et. al [56] also provide another point of 
comparison for the present numerical results. For their study, they looked at the effect of 
injector Mach number, pressure ratio, and nozzle divergence angle on jet primary wavelength 
and the Mach disk width and standoff distance. Jets were tested at Mach numbers ranging 
from 1 to 3 and for each case of Mach number the pressure ratio was varied over a wide range. 
Not only did they consider the characteristics of jets issuing into a quiescent medium, but 
they also considered cases of jets issuing into a supersonic coflow. Only the quiescent results 
are considered here, however. The quality and detail of the results in these experiments 
provide another good comparison point for the present numerical results. 

Converged results from VULCAN are compared to experimental results in Figures 15 
to 17. Each case consisted of a sonic, 300 K, 100% hydrogen jet issuing into quiescent air 
at an ambient temperature of 300 K and at a pressure of 2000 Pa. Jet pressure ratios of 
2, 3.8, 5, 10, 16.8, 20, and 27.6 were simulated. The non-integer pressure ratios correspond 
to experimental conditions. Hydrogen was chosen as the injectant due to the desire to use 
the same injectant gas as will be used in the later numerical simulations. This discrepancy 
in injectant gas between simulation and experiment is not anticipated to cause an issue, 
however, because the important parameters for this problem as identified by previous re- 
searchers are jet Mach number, pressure ratio, and ratio of specific heats [25, 57, 58]. For 
example, a semi-emperical relation by Billig et al. [25] expresses the Mach disk standoff 
distance in a quiescent crossflow as a function of jet Mach number and pressure ratio as 


In the same paper, Billig validates this relation against a wide range of experiments using 
different injectants. 



(1) 
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Figure 14: A plane bisecting the sonic jet is shown from the converged case of jet pressure 
ratio equal to 20. Parameters of the jet are shown, consisting of wavelength, w, Mach disk 
diameter, D m , exit diameter, D e , and Mach disk standoff distance x m . 
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Figure 15: Normalized Mach disk diameter as a function of jet pressure ratio. Numerical 
results from Vulcan are compared against the results of Wilkes et. al [55] and Love et. al 

[56]- 

Figure 15 compares computations to experiments of normalized Mach disk diameter 
over pressures ranging from roughly 5 to 40. The scatter of the results of Wilkes et. al can 
be explained due to the variation in Reynolds number for different runs. Regardless, the 
semi-empirical results of Love et. al at the same Mach number pass through the center of 
the NO-PLIF data, indicating general agreement between the two experiments. The present 
numerical results are indicated by solid triangles and show that the present results from 
VULCAN lie in the experimental scatter. 

Figure 16 shows the experimental and empirical jet Mach disk standoff location com- 
pared to the predictions by VULCAN for a sonic jet over the same range of pressure ratios. 
Here some divergence can be seen between the two experimental results as pressure ratio is 
increased. Love [56] does mention that the semi-empirical relation tends to overpredict the 
wavelength for high jet pressure ratios at low Mach numbers; however, the relation derived 
by Love is fit for the sonic case and matches their experimental data closely. Therefore, the 
cause of the discrepancy is not clear. The semi-empirical relation of Billig from Equation 
1 is also plotted and shows good agreement with the present results and the data taken by 
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Figure 16: Normalized Mach disk standoff distance as a function of jet pressure ratio. 
Numerical results from Vulcan are compared against the results of Wilkes et. al [55] and 
Love et. al [56]. 

Wilkes. Billig [25] also noted that Love’s relation for Mach disk standoff location tended to 
not match a general data set of experimental results. 

Figure 17 compares the experimental jet wavelength to that predicted by VULCAN 
at the same range of pressure ratios. As with the Mach disk standoff distance, there is 
a divergence between the two experimental results. As previously stated, the cause of 
this discrepancy is unclear. The present numerical results, however, agree well with both 
experimental results, lying between them. The results also lie within the experimental 
scatter. 

It has been shown that VULCAN is able to successfully predict the physics of sonic jets 
issuing into a quiescent medium. Numerical results have been compared against experi- 
mental data and semi-empirical relations for the parameters of jet wavelength, Mach disk 
standoff distance, and Mach disk width. The present numerical results accurately predict 
the experimental results from the literature. 
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Figure 17: Normalized jet wavelength as a function of jet pressure ratio. Numerical results 
from Vulcan are compared against the results of Wilkes et. al [55] and Love et. al [56]. 

2.2.2 Transverse Injection into a Supersonic Crossflow 

The ability of VULCAN to predict the physics of sonic injection into a supersonic crossflow 
will now be summarized. As with with injection into a quiescent medium, this case has 
been studied extensively by previous researchers at low Mach numbers. Therefore, this is 
another good case to study in order to validate VULCAN’s predictive capability in the 
flow regime of interest. As with the previous study, validation is occurring by injecting a 
sonic jet transverse to the flow. Although the cross flow Mach number of 2.1 is much lower 
than the cross flow Mach number in this dissertation, it is assumed that the results of this 
validation study can be extrapolated to higher Mach numbers due to the momentum-based 
nature of the jet trajectory [59] (i.e. no emergent effects due to increasing Mach number 
arises compared to an already supersonic cross flow). 

As with the injection of a jet into quiescent flow, the location of the Mach disk for a 
supersonic crossflow has typically been a parameter of interest. Cohen et al. [58] proposed 
the relationship for Mach disk midpoint height as a function of jet and freestream parameters 
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(2) 


D e Jj Mj [p/j + 1) [(l- 7oo) + 2700 M% 0 _ 


which was validated against experiments at conditions ranging from 0.55 < M a < 2.92 and 


1.0 < Mj < 2.2. 


The trajectory of the jet as it travels downstream has also been of interest to researchers. 


The penetration of the fuel is defined to be either the maximum point away from the wall 
where fuel is detected or the centerline of the jet. Schetz and Billig [57] , who chose the 
latter definition, used the curve fit 


to describe the trajectory of the jet centerline at different dynamic pressure ratios. This 


recently by Portz and Segal [60]. 

The grid topology used in this validation study was the same used in the study of a 
sonic jet issuing into a quiescent medium. While the pressure ratio was the main parameter 
of interest in the previous study, the current parameter of interest is the dynamic pressure 
ratio, J, of the fuel jet to the main stream. The dynamic pressure ratio was variable in this 
study while all other parameters were held constant. The Mach 2.1 cross stream is 300 K, 
2000 Pa air for all cases. The sonic fuel jet is 300 K hydrogen with a variable density to 
match the desired dynamic pressure ratio. 

The trajectory traveled by the fuel jet will now be examined. Figures 18 and 19 show 
the numerically solved hydrogen mass fraction on the flow centerline. Superimposed on the 
VULCAN solutions are the analytical relations of Equation 3 with different values of m and 
n. It is observed that the jet plumes simulated using VULCAN tended to underpredict the 
Schetz jet penetration, but in general followed the empirically-derived trend. 

The simulated location of the Mach disk will now be investigated and compared to 
Equation 2 using the same jet and crossflow parameters as those previously stated. Also 
included are experimental results by Schetz at similar freestream and jet Mach numbers. 
The VULCAN results are plotted and encompass a range of dynamic pressures from 4.75 



(3) 


relation has been compared successfully against experiment and has been validated more 
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Figure 18: The centerline hydrogen mass fraction field of a converged VULCAN case 
where J = 4.75 is shown compared to the analytical relation from Equation 3 with m = 1, 
n = 2.55. 



Figure 19: The centerline hydrogen mass fraction held of a converged VULCAN case where 
J = 16.35 is shown compared to the analytical relation from Equation 3 with m = 1.3, 
n = 3. 
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Figure 20: Vertical displacement of the Mach disk midpoint as a function of dynamic 
pressure. 

to 50. It is observed that in general VULCAN overpredicted the location of the Mach disk 
midpoint, although the error is limited to between 10-20%. 

It has been shown that VULCAN is able to successfully predict the penetration char- 
acteristics of sonic jets issuing into a supersonic crossflow. Numerical results have been 
compared against both experimental results and empirically- validated relationships. Over- 
and under-prediction was evident in the comparisons, but the errors are deemed to be ac- 
ceptable at the current level of validation fidelity; however, future attention should be given 
toward reducing these discrepancies. 

2.2.3 Axial injection into a Supersonic Coflow 

The ability of VULCAN to predict the physics of axial injection of fuel into a supersonic 
coflow will now be summarized. Numerical validation will be achieved by simulating a 
fuel jet issuing from the base of an expansion ramp situated in a supersonic airstream. 
Successful reproduction of experimental results for this case will show that VULCAN is 
capable of predicting the physics salient to such flows, including the penetration of the 
fuel jet, vorticity produced by the injector ramp, and the diffusion of the plume into the 
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Figure 21: A schematic of the Waitz expansion ramp injector is shown (from [62]). 

surrounding air. The injection and crossflow Mach numbers of this validation study are 
similar to those of ramp and strut injectors in the primary study, and therefore is intended 
to be a more direct validation compared to the previous two validation studies. 

The geometry used in this validation study was the same as evaluated by Waitz, et al. 
[61] in his study of Mach 1.7 helium injected into Mach 6 air from the base of an expansion 
ramp. The expansion ramp geometry is shown in Figure 23, the salient features being: a) 
the injector itself which is a wedge that turns the flow and injects fuel from its base, and 
b) a wall that turns away from the flow on either side of the ramp in order to create a 
shock wave at the ramp base to enhance mixing due to baroclinic torque between the fuel 
plume and the surrounding air. Vorticity is also created in the airstream as the fluid flows 
off either side of the ramp. The ramp studied here is unswept (i.e. parallel sides), although 
enhanced vorticity can be obtained by a swept ramp concept (i.e. sides converging toward 
the base). 

For the validation study, the geometry was set up to match the Waitz geometry, the 
experimental apparatus of which is shown in Figures 21 and 22. The injector had a height 
at its base of 25.4 mm, a width of 6.35 mm, and a length of 152 mm. The top of the 
ramp injector had a slope of 4.76° and the expansion ramp was set to the same angle. The 
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Figure 22: A photograph of the Waitz expansion ramp injector experimental apparatus is 
shown (from [62]). 

incoming Mach 6 air had a total pressure of 6.0 MPa and a total temperature of 520 K. 
The Mach 1.7 helium injectant had a total pressure of 23 kPa and a total temperature of 
223 K. These conditions corresponded to the baseline case conditions given by Waitz. 

A view of the computational geometry used in this study is shown in Figure 23. Symme- 
try planes were defined on the injector centerline and on the midpoint between neighboring 
injectors. All walls were defined to have a constant temperature of 300 K. During the 
course of his experiments, Waitz tripped the incoming boundary layer with gas so that it 
would have a controllable thickness. The baseline boundary layer thickness was 20% of the 
injector height. In order to provide the same thickness for the numerical study, a separate 
simulation was conducted on a 2D flat plate at the Mach 6 air conditions in order to allow 
a boundary layer to grow. At the plane where the boundary layer was 5.08 mm thick, the 
plane was saved in order to be extruded across the inflow of the primary simulation. 

A view of the computational grid in the vicinity of the injector is shown in Figure 
24. While greater grid density can be achieved for reasonable simulation times, it was 
decided to emulate the general grid layout and point density chosen by Waitz. Waitz used 
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Figure 23: The expansion ramp geometry of Waitz for the validation study is shown. 



Figure 24: The computational grid for the validation study in the vicinity of the injector 
is shown. 


point blanking 3 when building the injector geometry, something not employed in this study. 
Instead, the expansion ramp was designed to begin with a pole domain. This leads to a 
total of 334,000 grid points in the current topology instead of 266,000 in Waitz’s grid. 

The selection of turbulent Schmidt number has been shown to have a significant effect 
on the behavior of high speed mixing systems. The turbulent Schmidt number is a ratio of 
the turbulent viscous diffusivity of the system to the turbulent mass diffusion, defined by, 


c Vt 

SC ‘=D, 


(4) 


Typically lower turbulent Schmidt numbers will enhance the turbulent mixing present in 
the system while higher numbers will produce a more conservative estimate of the mixing 


3 Point blanking is a method of constructing a grid topology where instead of building the grid around 
the geometrical features, the grid points that occupy the boundary of the geometry are simply removed. 
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present in a system. Because this validation study endeavors to compare the simulated 
helium jet plume characteristics to experiment, it is evident that the Schmidt number is 
a key parameter. Therefore, a range of values were chosen between 0.1 and 1.0 and held 
constant through the entire domain. Recent studies have investigated the use of variable 
Schmidt numbers [63], but that is beyond the scope of this research. 

A unique solution was computed for each value of turbulent Schmidt number with 
a standard Menter k -uj turbulence model and compared against experiment to see what 
turbulence parameters produce the best results.. This process is known as calibration, 
which AIAA [64] defines as 

“Calibration: The process of adjusting numerical or physical modeling param- 
eters in the computational model for the purpose of improving agreement with 
experimental data” 

Once each solution reached convergence, cross-section planes were interrogated for the maxi- 
mum helium mass fraction at locations of 0, 1, 4, 8, and 13 x/hj (= x) distances downstream 
of the injector exit, where hj is the height of the injector base. These values were chosen 
because the experiment of Waitz readily had data available at these locations. The com- 
parison of Waitz’s experiment to each solution is shown in Figure 25. It is observed that 
simulations with lower turbulent Schmidt number had increased jet plume diffusivity and 
the scatter of the results indicates that the results of the simulation were highly sensitive 
to the chosen value of turbulent Schmidt number. Of all the values simulated, Set = 0.25 
appears to have the best agreement with experiment, overpredicting maximum helium mass 
fraction at x = 8 by 16% and underpredicting by 33% at x = 13. 

Results from the work of Waitz in Figure 26 show experimental and numerical helium 
mass fractions at x = 4, 8 and 13. It should be noted that the numerical simulations 
conducted by Waitz did not use any turbulence modeling. Nonetheless, the general shape of 
the plume at each location shows agreement between experiment and numerical simulation. 
The present numerical results carried out in VULCAN are shown in Figure 27 for x = 4, 
8, and 13 and Set = 0.25, 0.5, and 1.0. The helium contours plotted are the same values 
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Figure 25: Comparisons of maximum helium mass fraction as a function of normalized 
distance from the injector exit between experiment and numerical simulations at various 
turbulent Schmidt numbers. Lines are for visual aid only and do not imply a curve shape. 

as plotted in the experimental results of Waitz in Figure 26. Again, it is observed that 
good morphological agreement between the numerical and experimental cross sections were 
achieved. It is also evident that for each cross-sectional plane, the chosen value of Set had a 
weak effect on the shape of the plume; rather, the effect is most evident in the distribution 
of helium in the plume. Lower values of Set produced shallower gradients of helium mass 
fraction due to the greater turbulent molecular diffusivity of those cases. Conversely, high 
values of Set caused higher gradients of helium mass fraction due to inhibited diffusivity. 
It should also be noted that Set = 0.25 produced helium mass fraction contours that are 
similar to the experiment by Waitz, thus giving credibility for the choice of Set for this set 
of experimental conditions. 

It has been shown that VULCAN is able to successfully predict the physics of supersonic 
axial injection into a supersonic coflow. Many different mechanisms are present in such a 
flow, for example vorticity generated due to the ramp injector, baroclinic torque due to the 
passage of a shock wave through a helium jet surrounded by air, and diffusion and mixing 
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Figure 26: Experimental contours of helium mass fraction at normalized x coordinates of 
(a) 4, (b) 8, and (c) 13 from the injector exit are shown compared to numerical contours at 
normalized x coordinates of (d) 4, (e) 8, and (f) 13. 
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(g) (h) (i) 


Figure 27: Numerically solved contours for helium mass fraction are shown for (a) - (c) 
Set = 0.25, (d) - (f) Set = 0.5, and (g) - (i) Set = 1.0. The columns sequentially show 
normalized x coordinates of 4, 8, and 13 from the injector exit. 
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of helium into air. The present numerical results were compared against past experiment 
results and showed reasonable agreement between the two. Furthermore, the significant 
effect of turbulent Schmidt number was observed for this study, although the best value of 
Set for these conditions does not necessarily hold true for disparate conditions. 

2.3 Modeling and Simulation 

2.3.1 Performance 

2. 3. 1.1 Mixing Efficiency 


Mixing efficiency was used in this study as a measure of the effectiveness of a particular 
fuel injection concept to provide uniform mixing. Ideal conditions for combustion consist of 
uniformly distributed fuel- air at a local equivalence ratio of one at each point. Here mixing 
efficiency is an integrated value over a given flow plane and is defined by the ratio of mixed 
oxidizer mass flow rate divided by the total fuel plus oxidizer flow rate. The definition given 
by Mao [ 65 ] is 

mc>2,mix Ja,y=o Y/ipudA 

3 m = = — r 77 — JT ( 5 ) 

^02, total Ja,y= o ^ pudA 


where 


Yr={ 


Y for Y < Y s 


( 6 ) 


T EfY s for Y > Y s 

where Yg is the mass fraction of oxygen at stoichiometric proportions with hydrogen and 
air. Because the forebody mixture was globally fuel rich due to spillage of predominantly 
oxygen-containing flow from the open geometry, the mass fractions in these equations are 
that of oxygen; otherwise fuel mass fractions would be used. Mixing efficiency represents an 
upper bound on the combustion efficiency, which is a measure of fuel consumption rather 
than mixing. 


2. 3. 1.2 Stream Thrust Potential 


The typical loss figure of merit used to evaluate forebody flows is total pressure recovery. 
This is due to the ease with which it can be calculated or one-dimensionalized from a 
multidimensional flow. The use of such a parameter as a true indicator of irreversibility 
generation in the forebody flow implicitly assumes a nearly constant total temperature. For 
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forebody mixing problems, however, such an assumption breaks down due to the mixing of 
two streams of distinct total temperature and because flow spillage may be present. This 
issue was remedied by using the stream thrust potential parameter. 

Past studies have used thrust-based methods to identify, account for, and optimize on 
different loss mechanisms occurring in a flowpath [40, 66, 67, 68]. Its usefulness has been 
demonstrated in the design of complete engine flowpaths [69]. Although the current study 
only considers a single component of a full flowpath (i.e. a forebody injector), stream thrust 
potential still provides a good metric for comparing different injection concepts. The stream 
thrust potential is given by the equation 

F — TTlU e xp T PexpA-exp (7) 

where u exp and p exp are the velocity and static pressure when isentropically expanded from 
the simulated cross-sectional area to the expanded cross-sectional area A exp . The stream 
thrust potential is obtained by one-dimensionalizing the flow at a given plane according 
to a method that conserves mass, energy, and entropy [70]. Once the one-dimensionalized 
parameters are obtained, the flow at a given flow station is expanded isentropically to a 
reference nozzle area, shown schematically in Figure 28 for channel flow. For this study the 
reference nozzle area was set equal to the capture area of the notional flight vehicle. Mixing, 
heat loss, friction, and shock waves in general reduce stream thrust potential while heat 
addition and axial fuel injection increase stream thrust potential. Hence, as with mixing 
efficiency, using stream thrust potential by itself as a metric for comparison is insufficient 
because enhanced mixing will negatively impact stream thrust potential. 

2.3.2 Flow Conditions 

The simulated flow was on a flat plate at conditions representative of the forebody of an 
airbreathing vehicle on a typical constant dynamic pressure (71,820 Pa) trajectory at Mach 
12. Using the 1976 US standard atmosphere, such a point on the flight trajectory occurs at 
an altitude of 33.33 km. Representative flat plate conditions were obtained by processing 
the freestream flow through a five degree turn followed by a six degree turn. The notional 
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Figure 28: A schematic of the stream thrust potential methodology. 



Figure 29: The baseline forebody geometry used in this study is shown. 


forebody geometry used in this study to obtain the flat plate conditions is shown in Figure 
29. Freestream and flat plate conditions are summarized in Table 2. 

Each injector was sized assuming the engine flowpath captures 4.5 kg/s of air (corre- 
sponding to a small-scale engine) at the conditions in Table 2 and that the injector was one 
of five on the forebody combining to provide a globally stoichiometric fuel flow rate. For 
the case of transverse fuel injection, sonic injection with a dynamic pressure ratio of 1.0 
was ensured. For transverse injection, the fuel jet was underexpanded in order to prevent 
oblique shock structures from forming within the jet [71]. For ramp and strut injection, 
pressure matching between the injector exit and the inflow was assumed in an attempt to 
preserve vortical flow structure. An additional ramp injection case with sonic injection was 
simulated to evaluate the effect of higher pressure on injection and mixing. In all cases it 
was assumed that the fuel system is able to supply the pressures required for the injection 
system. The assumed and derived injector conditions are summarized in Table 3. 

The inflow boundary layer thickness was set to be 10 mm in this study. To accomplish 
this, a separate flat plate simulation was conducted at the conditions in Table 2. At the 
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Table 2: Freestream and flat plate conditions. Flat plate conditions are obtained by shock- 
ing the freestream flow through a 5 degree turning followed by a 6 degree turning. All cases 
use an air molar composition of 79% N 2 and 21% O 2 



Freestream 

Flat Plate 

Static pressure 

713 Pa 

9111 Pa 

Static temperature 

232 K 

541 K 

Mach number 

12.0 

7.67 


Table 3: Fuel injector exit conditions. All cases use 100% H 2 at a static temperature of 
390 iy, 



Transverse 

Ramp 

(Matched) 

Ramp 

(Sonic) 

Strut 

(Matched) 

Injector Mach number 

1.0 

3.92 

1.0 

3.92 

Dynamic pressure ratio 

1.0 

0.26 

0.07 

0.26 

Static pressure 

535,000 Pa 

9111 Pa 

35,715 Pa 

9111 Pa 

Injector exit area 

52.6 mm 2 

790 mm 2 

790 mm 2 

790 mm 2 


location where the boundary layer thickness grew to 10 mm, the flow profile was output at 
that point and extruded across the inflow plane of the injector simulations. 

2.3.3 Injector Geometries 

2.3.3. 1 Flush- Wall Injection 

In this study transverse injection was simulated through flush-wall nozzles of varying exit 
geometry. While flush wall injectors could be angled to achieve an axial component of the 
fuel flow, only transverse injection was considered herein in order to understand the limits 


Diamond Wedge 



di 

Figure 30: Notional geometry for flush wall injectors. 
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Table 4: Injector geometry for flush-wall injection. All injectors have the same exit area. 


Geometry 

di (mm) 

d 2 (mm) 

Circular 

8.19 

- 

Diamond 

24.44 

4.31 

Wedge 

17.28 

6.09 

Stinger 

36.83 

1.64 


of mixing obtainable by a flush wall injector. Injector concepts to be considered included a 
circular wall injector as well as injectors with leading edge geometrical modifications that 
make them more suitable for a forebody injection environment. 

The flush wall nozzle geometries studied herein are depicted in Figure 30 with parameters 
given in Table 4. Of the four geometries considered, three are primitive shapes designed to 
give a basic understanding of the effect of leading- and trailing-edge bluntness and sharpness 
on the evolution of mixing and losses downstream of the injection location. Sharp leading 
edge modifications to flush-wall injectors were proposed to reduce upstream recirculation of 
fuel. The fourth geometry was proposed by Hirano [30] and was demonstrated to produce 
enhanced plume area and penetration in a low Mach number cross flow. However, it was 
not characterized in terms of mixing efficiency or losses. This geometry, referred to as a 
stinger, was evaluated to study its effect in high Mach number cross flow. 

2. 3. 3. 2 Ramp Injection 

The ramp geometry has been studied and used in internal supersonic flows for the axial 
injection of fuel as well as inducing vortical flow structures to enhance mixing. Marble [34] 
proposed an expansion ramp concept in order to enhance mixing through the interaction of 
the fuel plume with the shock wave emanating from the base of the ramp. This interaction 
induces baroclinic torque due to the unaligned density and pressure gradients at the fuel-air 
boundary and causes a roll-up of the fuel in the direction the shock sweeps through the fuel 
plume. Ramp injectors that are swept induce additional vorticity as flow rolls off of the 
edges of the ramp [26] . 

In this study the expansion ramp injector of Waitz [72] was modified for the specific 
requirements of forebody fuel injection. The face of the injector base below the injection 
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Side View 



Figure 31: Side-, top-, and perspective-view geometry for modified ramp injectors. 


Table 5: Modified ramp cases. 


Case 

P (deg) 

Base Wedge? 

Swept ramp (no wedge) 

3.0 

No 

Swept ramp (with wedge) 

3.0 

Yes 

Unswept ramp (no wedge) 

0.0 

No 

Unswept ramp (with wedge) 

0.0 

Yes 


location was cut back upstream in order to eliminate the fueled recirculation region that 
would otherwise be observed at the base. Such a recirculation region increases the possibility 
of autoignition on the forebody due to the hot stagnation region behind the ramp. The 
geometry modification was done in two ways. The first was by simply cutting a cavity out 
of the base of the injector below the duct that is supplying fuel. This type of geometry has 
been studied by Parent et al. [35, 24, 36, 23, 37, 38] at similar conditions and is referred to 
as a cantilevered ramp injector. The second modification cut back the injector base in an 
expansion wedge, as shown in Figure 31. This geometry was used in an attempt to increase 
pressurization on the base of the injector while discouraging recirculation. A summary of 
the cases tested are given in Table 5. 

2. 3. 3. 3 Strut Injection 

Like the ramp injector, the strut has also been used in internal flows for supersonic fueling 
and mixing. The typical strut extends through the entire flowpath height and is character- 
ized as a fuel placement device rather than a fuel mixing device. The strut injector may 
potentially be the ideal fuel injector for large-scale (e.g. orbital access) flowpaths where the 
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Side View 



Figure 32: Side- and top-view geometry for strut injector. 


ability for ramps or flush-wall injectors to access the core flow is diminished. Geometries 
that have heights deep into the fuel core but do not access the entire flowpath height have 
been called either struts or pylons in the literature — for the purposes of this study they will 
also be called struts. Up to this point studies of forebody injection via struts have not been 
deeply investigated in the literature. 

This study considered a single strut injector geometry that was swept back 45 degrees 
and injected fuel away from the injector centerline at a 10 degree angle. A schematic of the 
injector is shown in Figure 32. As with the ramp injector, the strut injector was placed on 
an expansion ramp in order to take advantage of shock-induced baroclinic torque between 
the fuel plume and the surrounding air. Rectangular ports supplied fuel injection with 
expansion wedges cut between injectors in order to reduce base recirculation and increase 
pressurization. The strut design proposed here is a non-optimized point design meant to 
be a pathfinding investigation into forebody strut injection. 

2.3.4 Computational Domain 

The computational grids had approximately 5.8 - 7.8 million grid points with clustering 
near the injector location. The computational domains were representative of one-half of 
the physical domain by using a symmetry boundary condition along the flow centerline 
bisecting the injector. All walls were set to be isothermal surfaces at a temperature of 
500 K. The inflow boundary was fixed at the conditions given in Table 2, representative of 
conditions on the second ramp of the notional Mach 12 vehicle forebody. Shocks that would 
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Figure 33: A side- and top-view schematic of the computational domain is shown. 



Figure 34: The circular injector is shown positioned in the portion of the computational 
domain where performance sampling occurs. 


otherwise appear due to forebody flow turning were ignored for this study. A 10 mm portion 
of each fuel injection tube was simulated with slip walls; thus, the effects of boundary layer 
buildup and vorticity production within the injector nozzle were neglected. A schematic 
of a general computational domain with a flush wall injector is shown in Figure 33. The 
dimensions of the domain are representative of the second ramp of the Mach 12 reference 
vehicle. In each case there was a 50 mm ramp leading up to each injector, which occupies 
a 713 mm long region until the notional cowl closure location. The distance between the 
injector centerline and the halfway point between neighboring injectors was 34.5 mm. 

Structured grids were used throughout each topology due to VULCAN only accepting 
such grids at the time of this study. The topologies upstream and downstream of the 
injectors were simple H-grids with spacing at the wall in each case designed to be 0.127 
mm. Grid spacing normal to the wall was designed to gradually grow from 0.07% at the 
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Wedge Injector 



Figure 35: The wedge injector is shown positioned in the portion of the computational 
domain where performance sampling occurs. 



Figure 36: The diamond injector is shown positioned in the portion of the computational 
domain where performance sampling occurs. 
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Figure 37: The stinger injector is shown positioned in the portion of the computational 
domain where performance sampling occurs. 


49 



Figure 38: The unswept ramp is shown positioned in the portion of the computational 
domain where performance sampling occurs. 



Figure 39: The strut injector is shown positioned in the portion of the computational 
domain where performance sampling occurs. 
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Figure 40: The computational grid in the vicinity of the circular injector is shown. 



Figure 41: The computational grid in the vicinity of the wedge injector is shown. 


wall to 2.75% away from the wall according to the hyperbolic tangent stretching function. 
These characteristics were chosen to target a y + near 10 upstream of each injector. This 
value of y + was deemed appropriate due to the use of wall functions in VULCAN. 

The grid topologies in the vicinity of the flush- wall injectors are shown in Figures 40 to 
43. The circular injector consisted of a central H-grid wrapped with a C-grid that was then 
connected to the rest of the wall topology consisting of H-grids. Thirty grid points were 
placed along the diameter of the circular injector with the resulting injector grid spacing 
used to place grid points in the diamond, wedge, and stinger geometries. For the injectors 
with sharp leading or trailing edges, H-grids with pole domains were used to collapse the 
grid into single point (s) due to the shallow angles involved. 

The grid topology near the unswept ramp with an expansion wedge feature is shown 
in Figure 44 and is representative of the topology of the four ramp geometries considered. 
Pole domains were used to run the H-grids on the side of the injector toward the single 
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Figure 42: The computational grid in the vicinity of the diamond injector is shown. 



Figure 43: The computational grid in the vicinity of the stinger injector is shown. 


52 









Figure 44: The computational grid in the vicinity of the unswept ramp wedge injector is 
shown. 



Figure 45: The computational grid in the vicinity of the strut injector is shown. 

point on the leading edge of the injector. For cases with the expansion wedge present, an 
H-grid with pole domain was similarly used on both the wall and underside of the injector. 
The injector exit of the ramp was designed to have 30 grid points across both the width 
and the height of the exit. 

The grid topology near the strut injector is shown in Figure 45. The topology of the 
strut injector used H-grids with pole domains being located at the leading edge of the strut 
as well as on the expansion wedge between the two injection ports. 

A summary of the total number of grid points in each computational domain and in the 
performance sampling region shown in Figure 33 are given in Tables 6 and 7, respectively. 
Downstream of each injector the I, J, and K directions were aligned with the X, Y, and Z 
coordinates, respectively. In the performance sampling region each topology was sized to 
give roughly the same number of grid points in the wall-normal and lateral directions while 
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Table 6: Number of grid points full computational domain for each injector concept. 


Total Cells 

Circular 

6.6M 

Diamond 

6.2M 

Stinger 

6.4M 

Wedge 

5.8M 

Ramp 

7.8M 

Strut 

5.9M 


Table 7: Number of grid points in the performance sampling region in each direction as 
well as total cells. 



I 

J 

K 

Total Cells 

Circular 

541 

69 

121 

4.4M 

Diamond 

531 

60 

121 

3.8M 

Stinger 

522 

60 

121 

3.7M 

Wedge 

537 

60 

121 

3.8M 

Ramp 

396 

73 

124 

3.5M 

Strut 

434 

57 

121 

2.9M 


the longitudinal spacing was varied depending on the downstream distance occupied by the 
injector itself. 

2-4 Mixing Results 
2.4.1 Flush Wall Injection 

This subsection investigates the relative performance between different transverse wall in- 
jectors at the same fuel injection mass flow rate and injector area. Figure 50 shows fence 
plots of hydrogen mass fraction contours for the flush wall injectors. In general, the fuel 
plumes from the flush wall injectors remained close to the wall. This is undesirable when 
trying to prevent autoignition of the fuel in the boundary layer. Even the stinger, which 
exhibited the best penetration of fuel into the main flow of the four cases, maintained a fuel 
rich region near the wall. The circular injector had a large recirculation region upstream 
of the jet which increases the possibility of autoignition due to the increased residence time 
of entrained fuel combined with high temperatures in this region. Switching to a geometry 
with a sharp leading edge largely eliminated this upstream recirculation region, which is 
beneficial for mitigating autoignition. 

Contours of wall y + values for each flush wall injection geometry are shown in Figures 
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Figure 46: Contour plots of wall y + for the (a) full flowpath and (b) in the vicinity of the 
fuel injector. Flow is from the left to the right and hydrogen issues from the injector out of 
the page. 


Table 8: Maximum y + values along the wall for each flush-wall fuel injection concept. 


Circular 

Diamond 

Stinger 

Wedge 

Vtnax 20.13 

63.93 

54.50 

55.19 


46 to 49 using the same color scale with the maximum values for each case given in Table 8. 
Upstream of each injector a y + of 11 was achieved, while downstream in the mixing region 
values dropped below 10. The maximum values for each case occurred in the vicinity of the 
injectors near the bow shock created by the fuel plume. Although these maxima occur in 
small regions relative to the full flowpath, they are still within the valid range of using wall 
matching functions. 

Figure 51 shows the performance of each injector as a function of downstream distance 
for the metrics of mixing efficiency and stream thrust potential. Note that the spatial 
coordinate has its origin at the leading edge of each injector. The mixing efficiency curve 
for each injector followed a monotonically increasing trend as shown in Figure 51(a). In 
terms of mixing efficiency, the stinger injector was the best performer by a significant margin 
while the classic circular injector and the diamond injector were the worst performers at 
cowl closure. It is also notable that the stinger injector has the greatest fineness ratio versus 
the circular injector. Because each injector type had only one example of fineness ratio (i.e. 
injectors with the same shape but different aspect ratios were not studied), the ability to 
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Figure 47: Contour plots of wall y + for the (a) full flowpath and (b) in the vicinity of 
the diamond fuel injector. Flow is from the left to the right and hydrogen issues from the 
injector out of the page. 



Figure 48: Contour plots of wall y + for the (a) full flowpath and (b) in the vicinity of 
the diamond fuel injector. Flow is from the left to the right and hydrogen issues from the 
injector out of the page. 
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Figure 49: Contour plots of wall y + for the (a) full flowpath and (b) in the vicinity of 
the diamond fuel injector. Flow is from the left to the right and hydrogen issues from the 
injector out of the page. 


identify this as the dominant driver of mixing efficiency was limited. 

Stream thrust potential as a function of distance from the injector leading edge is shown 
in Figure 51(b). Each of the injectors with sharp leading edges had roughly the same stream 
thrust potential loss with the stinger injector being the worst of the three as the cowl 
closure location was approached. For most of the mixing region the circular injector had 
a worse stream thrust potential compared to the sharp leading edged injectors and better 
performance as the cowl closure location was approached. In general, the stream thrust 
potential experienced a slight decline initially due to friction and mixing. At around x = 
0.2 m the stream thrust potential experienced a sharp decline due to flow spillage from the 
top of the investigated domain. Because the mass flow rate is a dominant parameter in the 
stream thrust potential equation this caused a corresponding drop in stream thrust potential 
(even if the studied area and total pressure remained otherwise unaffected). Therefore flow 
spillage was accounted for in the stream thrust potential calculation (as opposed to total 
pressure loss, which may be unaffected by spillage). 

It is difficult to analyze mixing efficiency and stream thrust potential separately to get 
an idea of which injector performs the best. For example, just considering mixing efficiency, 
the stinger injector performed the best while stream thrust potential indicates that the 
circular injector was the best performer. Additional insight can be gained by considering 
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Figure 50: Fence plots of hydrogen mass fraction for the flush-wall injection cases. 


58 






(a) 



(b) 



(c) 


Figure 51: Performance of transverse injection concepts: (a) mixing efficiency versus x, (b) 
stream thrust potential versus x, and (c) stream thrust potential versus mixing efficiency. 
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Figure 52: Fuel equivalence ratio contours for the (a) circular, (b) diamond, (c) wedge, and 
(d) stinger injection case at x = 0.1 m and 0.713 m. Lines are plotted as <f> = 0.1, 1,2. ..7. 
Circular marker indicates vortex center. 


stream thrust potential as a function of mixing efficiency, plotted in Figure 51(c). This plot 
shows how much stream thrust potential was recovered for a given level of mixing for each 
concept. This curve gives the analyst an idea of the value of each injector. Among the flush 
wall injectors considered, the stinger injector was the best value because for any given level 
of mixing it had the greatest stream thrust potential. Because it was able to achieve higher 
levels of mixing overall, it would be able to mix out faster than any other concept (thereby 
reducing the required forebody length). 

Plots of fuel equivalence ratio just downstream of the injector and at the cowl closure 
location are shown in Figure 52. Just downstream of the injection location for each injector 
a vortex core formed in the fuel plume. The vortices at the cowl closure for each injector case 
were relatively weak. The vorticity observed in the stinger case combined with the rising 
fuel bubble visible near the centerline led to enhanced diffusivity of the fuel plume that 
yielded its advantageous mixing performance. Because of the importance of characterizing 
fuel spillage, the spillage values at cowl closure are given in Table 9. Low values of fuel 
spillage were observed for each injection concept. Because the displacement of the fuel 
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Table 9: Flush-wall injected fuel spilled at cowl closure as a percentage of the nominal fuel 
injection mass flow rate of 0.02639 kg/s. 


Circular 

Diamond 

Stinger 

Wedge 

% Fuel Spilled 2.9% 

1.4% 

3.0 % 

1.8 % 



y': 5 1015202530354045505560 



y': 5 1015202530354045505560 



(b) 

Figure 53: Contour plots of wall y + for the full flowpath of the unswept ramp fuel injector 
with cavity at both a (a) pressure matched and (b) sonic pressure condition. Flow is from 
the left to the right and the ramp itself is blanked out. 


plume is a function of dynamic pressure ratio, these low values of spillage may not hold for 
higher dynamic pressure ratios. However, it is shown that for a dynamic pressure ratio of 
unity the fuel is able to access much of the captured flow area while spilling only a small 
percentage. 

2.4.2 Ramp Injection 

This subsection assesses the performance of forebody ramp injection for the cases listed 
in Table 5. The objective was to assess the relative performance of ramp injectors as a 
consequence of varying ramp geometry and injector pressure. For each case the equivalence 
ratio was held constant and consideration was made toward the spillage of air from the 
capture streamtube. 

Contours of wall y + values for each simulated case are shown in Figures 53 through 56. 
The inflow y + value for each case was 11 with maximum values given in Table 10, each 
occurring downstream of the ramp in the mixing region and within the range of validity for 
using wall matching functions. 
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Figure 54: Contour plots of wall y + for the full flowpath of the unswept ramp fuel injector 
with expansion wedge at both a (a) pressure matched and (b) sonic pressure condition. 
Flow is from the left to the right and the ramp itself is blanked out. 
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(b) 

Figure 55: Contour plots of wall y + for the full flowpath of the swept ramp fuel injector 
with cavity at both a (a) pressure matched and (b) sonic pressure condition. Flow is from 
the left to the right and the ramp itself is blanked out. 


Table 10: Ramp injected maximum wall y + values. Legend: C(antilivered), W(edge), 
M(atched pressure), S(onic pressure) 


Base 

Unswep 

C 

t Ramp 
W 

Swept 

C 

Ramp 

W 

Pressure 

M 

S 

M 

S 

M 

S 

M 

S 

Vmax 

16.77 

20.37 

17.37 

20.99 

20.56 

15.05 

20.47 

14.88 
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(b) 

Figure 56: Contour plots of wall y + for the full flowpath of the swept ramp fuel injector 
with expansion wedge at both a (a) pressure matched and (b) sonic pressure condition. 
Flow is from the left to the right and the ramp itself is blanked out. 

Figure 57 shows fence plots of hydrogen mass fraction contours for the ramp injector 
without the base expansion wedge. The pressure matched ramps, both swept and unswept, 
were able to inject fuel such that it mixed away from the wall. Increasing the injection 
pressure to a sonic condition caused the fuel to expand toward the wall and, in the case 
of the swept ramp, created a recirculation region under the ramp that consequently mixed 
fuel near the wall. 

An initial evaluation of the performance of each injector was made using stream thrust 
potential and mixing efficiency. Figure 58 shows the effect of geometry on the mixing 
efficiency and stream thrust potential for a fuel plume pressure matched to the air inflow. 
Geometry considerations included the sweep of the ramp and the presence of a wedge on 
the underside of the fuel duct. Figure 58(a) shows the mixing efficiency as a function of 
x-location and begins at the injector exit. This figure demonstrates an advantage in mixing 
by the addition of the wedge feature and an advantage in mixing by ramp sweep. The latter 
result is expected due to the enhanced vorticity of the swept ramp case. 

A larger difference in performance is observed by considering Figure 58(b), which plots 
stream thrust potential as a function of distance for each ramp case. Here the unswept ramp 
had a consistent performance benefit over the swept ramp. This was due to the fact that 
the wider leading edge of the swept ramp caused overall increased mass flow spillage and 
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Figure 57: Fence plots of hydrogen mass fraction for the ramp injection cases without the 
base expansion wedge. 
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Figure 58: Performance of ramp injection concepts for matched pressure: (a) mixing effi- 
ciency versus x, (b) stream thrust potential versus x, and (c) stream thrust potential versus 
mixing efficiency. Dashed lines indicate cases with the expansion wedge feature. 
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Figure 59: Stream traces that indicate flow arcing. 
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Figure 60: Performance of ramp injection concepts for sonic pressure: (a) mixing efficiency 
versus x, (b) stream thrust potential versus x, and (c) stream thrust potential versus mixing 
efficiency. Dashed lines indicate cases with the expansion wedge feature. 
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Figure 61: Fuel equivalence ratio contours for the swept and unswept ramp cases at x = 
0.21 nr and 0.713 nr. Lines are plotted as <p = 0.1, 1, 2. ..7. Circular marker indicates vortex 
center. 
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flow expansion due to the shock wave originating from the ramp leading edge. Eventually 
the flow expansion turned the downstream velocity vectors near the cowl-closure y-location 
back toward the body and increased the mass flow rate in the studied region. This is 
demonstrated by plotting centerline streamtraces in Figure 59. This caused a corresponding 
increase in stream thrust potential that affected each case similarly (and therefore does not 
substantially affect the relative performance of each concept). 

Figure 58(c) plots the value of each injector concept by plotting the stream thrust 
potential as a function of mixing efficiency. Here it is observed that adding the expansion 
wedge provided a slight advantage over the cases without the feature. Furthermore, using an 
unswept ramp provided the best stream thrust potential for the amount of mixing obtained. 

Figure 60 plots performance metrics of the ramp cases when the injector pressure is 
raised by a factor of 3.92 (resulting in sonic injection). Figure 60(a) shows enhanced mixing 
on a case- by-case basis with increased fuel pressure. This is due to the underexpanded fuel 
plume reaching more of the main flow and inducing additional vorticity (discussed later). 
At cowl closure the high pressure, swept ramp case had 22% better mixing relative to the 
pressure matched case for the same geometry. However, the high pressure case had lower 
stream thrust potential relative to the pressure matched case, as seen in Figure 60(b). Be- 
cause of this, the stream thrust potential versus mixing efficiency curve in Figure 60(c) is 
lower than that for the pressure matched case. Therefore, among the cases considered, pres- 
sure matching is preferred-although this conclusion does not necessarily imply superiority 
over all other injector pressures. 

A comparison of cross sectional slices of the flow for the pressure matched and sonic 
ramps are shown in Figure 61. The cross sectional slices for the cases with an expansion 
wedge on the injector base are substationally the same as those without, so the wedge cases 
are omitted here. The pressure matched cases both had the classic “kidney bean” shape at 
the cowl closure location that is indicative of enhanced mixing due to vorticity. The sonic 
injector had an additional lobe on the periphery of the fuel plume due to a secondary vortex 
that further enhanced mixing over the pressure matched case. The pressure matched case 
began to form the second vortex at an intermediate location but it was quickly absorbed by 
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Table 11: Ramp injected fuel spilled at cowl closure as a percentage of the nominal fuel 
injection mass flow rate of 0.02639 kg/s. Legend: C(antilivered), W(edge), M(atched pres- 
sure), S(onic pressure) 
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the primary vortex. The sonic injection case also had enhanced mixing near the wall which 
may have consequences for the possibility of boundary layer autoignition. This infiltration 
near the wall is due to the fact that fuel was entrained upstream of the ramp injector exit 
for the higher pressure case and was mixed closer to the wall due to the vorticity initiated 
by the ramp itself. Fuel spillage was characterized for each ramp case and is summarized 
in Table 11. For each case the injected fuel only had small levels of spillage, thereby giving 
credibility to the ramp design in terms of being able to capture the injected fuel. 

2.4.3 Strut Injection 

Note that the discussion contained in this subsection for strut injection is only for a point 
strut design and does not necessarily represent an optimum configuration. In this disser- 
tation, parametric studies of strut sweep angle, fuel injection angle, and fuel pressure were 
not conducted for the given configuration. The purpose of this discussion is to get an idea 
of some of the general performance characteristics of a strut in forebody flow conditions as 
a starting point for future analysis. 

Contours of wall y + values are shown in Figure 62 for the strut injector. The inflow 
y + value was 11 with a maximum value of 50.11 occurring downstream of the strut on 
the centerline of the mixing region, within the range of validity for using wall matching 
functions. 

Figure 63 shows a fence plot of hydrogen mass fraction contours for the strut injector. 
Overall the fuel plume had excellent mixing into the surrounding air, although there was 
some undesirable mixing near the wall. This could be mitigated in future strut designs by 
placing fueling ports farther from the wall. The effect of shock waves turning fuel toward 
and away from the centerline can also be observed in this figure. 
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Figure 62: Contour plots of wall y + for the full flowpath of the strut injector. Flow is from 
the left to the right and the strut itself is blanked out. 
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Figure 63: Fence plot of hydrogen mass fraction for the strut injection case. 

Figure 64 summarizes the performance characteristics of the strut injector concept. 
Because the strut is a fuel placement apparatus, the flow was able to be seeded in multiple 
locations (four for this injector). The benefit to mixing efficiency is visible in Figure 64(a). 
The strut injector had significantly higher mixing efficiency at cowl closure than the ramp 
injector. One reason for this is due to the multiple-port characteristic of the apparatus. 
By having multiple ports, the contact surface between the fuel and the air was lengthened, 
allowing for improved mixing of fuel into air. 

The stream thrust potential of the strut in Figure 64(b) follows a monotonically decreas- 
ing trend. Near the cowl closure point there was a shock wave that spilled an additional 
amount of fuel. When considering the value of the mixing efficiency in Figure 64(c) it is 
evident that the strut was able to maintain a relatively high level of stream thrust potential 
for the mixing efficiency it achieves. The relative performance of this concept with respect 
to the other concepts is discussed in the next section. 

Cross sections of fuel equivalence ratio for the strut injector just downstream of the 
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Figure 64: Performance of the strut injection concept: (a) mixing efficiency versus x, (b) 
stream thrust potential versus x, and (c) stream thrust potential versus mixing efficiency. 
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Figure 65: Fuel equivalence ratio contours for the strut injector case at x = 0.16 and 0.713 
m m. Lines are plotted as <p = 0.1, 1, 2. ..7. Circular marker indicates vortex center. 

trailing edge and at cowl closure are shown in Figure 65. In general, the flow for the strut 
case was shock-dominated, with the fuel jet turned toward and away from the centerline as 
shock waves reflect off of the two symmetry boundary conditions. The fuel plume mixed 
well overall with a smaller core than found with previous injector concepts. There is also 
vorticity present near the bottom of the fuel plume that further enhanced mixing and is 
sustained through cowl closure. Fuel spillage for the strut injector at cowl closure was 3.0% 
relative to the nominal injected fuel mass flow rate of 0.02639 kg/s. This value of spillage, 
although low, may be reduced if the strut height is similarly reduced. 

2.4.4 Injector Comparisons 

Comparison between injection strategies was conducted by selecting the best value perform- 
ers from each injector category (transverse, ramp, and strut) and evaluating their relative 
performance. The best performers from the first two categories were the stinger and unswept 
ramp with expansion wedge. The strut is included for completeness. Figure 66 summarizes 
the performance of each top-performer. 

Figure 66(a) shows the mixing efficiency as a function of location for each injector con- 
cept. While the strut injector began to inject fuel farther downstream, it rapidly approached 
and exceeded the mixing efficiency achieved by the transverse stinger injector. The ramp 
injector had significantly lower mixing efficiency compared to the other two concepts, even 
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Figure 66: Performance of the best injection concepts: (a) mixing efficiency versus x, (b) 
stream thrust potential versus x, and (c) stream thrust potential versus mixing efficiency. 
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considering its delayed injection location. 

Figure 66(b) compares the stream thrust potential of each injector concept. The best 
comparison is made at the cowl closure location which represents the true captured stream 
tube for each concept. At this point, the ramp wedge maintained the most stream thrust 
potential with the transverse injector having performed the poorest. This is an expected 
result due to the axial direction of the ramp injection and the high loss nature of transverse 
injection. The strut injector performed worse than the ramp injector, which may be in part 
due to the enhanced mixing of the strut as they both had relatively similar thrust potential 
at their injection locations. 

Figure 66(c) plots the stream thrust potential as a function of mixing efficiency. The 
location of greatest mixing efficiency for each concept was at cowl closure. The ramp injector 
was found to perform between the strut and transverse injectors, with the strut having had 
the best stream thrust potential for the amount of mixing it achieved. The strut also implies 
a much shorter mixing length (similar to the transverse injector, only with much less stream 
thrust potential loss) over the ramp injector. This indicates the competitiveness of the strut 
injection concept and that it should be included in future forebody injection studies. 

2.4.5 Effect of Turbulent Schmidt Number 

The effect of turbulent Schmidt number on the performance metrics of the stinger, unswept 
ramp with wedge, and strut will now be discussed. It was observed in the validation study 
with the Waitz geometry that the selected value of Set had a strong effect on the diffusion 
between the injectant and the surrounding air stream. This will have an effect on the 
calculated mixing efficiency for each injector concept. The value Set = 1.0 was chosen with 
the underlying assumption that while the true value of the performance metrics will not 
be necessarily calculated, this was allowable because the primary interest was comparison 
between injector concepts. 

Three values of Set were compared- 1.0, 0.5, and 0.25 — with Set = 1.0 being the baseline 
case that has already been the basis of the comparisons described to this point. All other 
simulation parameters were held constant. Figure 67 shows the thrust potential plotted 
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as a function of mixing efficiency for the stinger, strut, and unswept ramp with wedge for 
each value of turbulent Schmidt number. It is observed that the selected value of Schmidt 
number does impact the level of mixing achieved for each concept, although this effect is 
varied depending on the injector studied. For example, the stinger injector maximum mixing 
levels varied between 0.42 and 0.46, compared to the strut that varying between 0.49 and 
0.58. This may be due to the enhanced interfacial area between the strut fuel plume and 
the surrounding air. Also of note is that the thrust potential levels were not substantially 
affected, as might have been expected due to the losses normally incurred due to mixing. 
For each case of turbulent Schmidt number, the trend and range of thrust potentials was 
roughly the same. This is due to the fact that the dominant contributor to loss was the 
spillage of flow from the domain. 

An important conclusion to draw from looking at the effect of turbulent Schmidt number 
is that when comparing different injector concepts, the same conclusions can be drawn as 
was drawn before. For example, the stinger performed the best when valued against mixing 
efficiency. It also achieved the greatest level of mixing and therefore will have the shortest 
mixing length of the considered concepts. Therefore, while the solution was sensitive to the 
selection of turbulent Schmidt number, the present comparative study was not substantially 
dependant on the selection of the parameter. 

It should be noted that this study cannot give any recommendation on which value 
of turbulent Schmidt number should be chosen for the present regime of flight based on 
the current results. This underscores the importance of producing experimental validation 
data in a regime relevant to the forebody of a hypervelocity vehicle. Such experimental 
databases will provide important information for calibrating future RANS simulations to 
produce results that can accurately predict losses and mixing achieved for forebody fueling 
concepts. 

2.5 Chapter Summary 

The current study of injection and mixing on a vehicle forebody at hypervelocity conditions 
concerns the efficiency of three major categories of injection: flush-wall, ramp, and strut 
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Figure 67: Comparisons of stream thrust potential versus mixing efficiency between the (a) 
stinger, (b) unswept ramp with wedge, and (c) strut injectors for Set = 0.25 - 1.0. 
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injection. Each injector was sized and designed to provide stoichiometric fuel to air mass flow 
given a notional forebody design and vehicle capture. In evaluating the value each injector 
provided to the system, two metrics were used: mixing efficiency and thrust potential. 
Mixing efficiency was used to gauge how well each concept mixed (and is a metric that is 
typically maximized) while thrust potential was used to gauge what losses were incurred 
through the injection and mixing process. It was proposed that comparing the two metrics 
directly gave a valuable means to compare different injection concepts. 

Flush-wall injection was shown to overall achieve high levels of mixing, but at the cost 
of very high losses. This was mostly due to the fact that a perpendicularly-directed fuel jet 
produced a bow shock ahead of the fuel plume that provided a large amount of disturbance 
to the main flow. Furthermore, because the fuel was issued normally to the upstream 
velocity vector, the system was unable to take advantage of the thrust produced by the 
jet. Worse still, energy must be expended to turn the fuel plume in the axial direction. 
While flush-wall, angled injection would be able to alleviate some of the losses incurred, it 
would be at the expense of mixing efficiency. In addition to the quantitative performance 
of the flush-wall injectors, other qualities of the flow proved undesirable. For example, a 
recirculation region tended to form upstream of the injector that entrained fuel and provided 
a long fuel residence time. This increased the risk of autoignition occurring in such regions, 
something that is needed to be avoided for forebody injection concepts. 
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CHAPTER III 


AUTOIGNITION MITIGATION 


The premixed or mixing flow on a forebody shares a cooling requirement with an internal 
flow in that surface temperatures and heat loads below the wall material limit are needed for 
basic structure survivability. However, such forebody flows have an additional requirement: 
to minimize (and hopefully eliminate) the presence of ignition until the entrance of the 
combustor is reached. Premature burning not only consumes fuel before the fuel/air mixture 
can reach the combustor, but it also introduces pressure drag on the vehicle overall due to 
elevated pressures on the forebody. Like conventional scramjets, shcramjets are likely to 
operate on a thin margin of thrust — the introduction of both pressure drag and thrust losses 
(both through premature fuel consumption) may create design and mission infeasibility. 

Both the film cooling and transpiration concepts have been studied a great deal in the 
literature and are used in many different aerospace systems encountering high heat load 
conditions. They have been considered and used for the protection of rocket nozzles [73], 
supersonic inlets [74], gas turbine blades [75, 76, 77], nose cones [78], reentry vehicles [79], 
and rocket combustors. They are generally distinguished in that film cooling works by 
injecting a cooled gas near the body through a tangential slot and transpiration cooling 
injects gas through perforations in the wall. Both concepts alter the thermal structure of 
the wall boundary layer as well as convect heat away from areas of high heat load. This 
introduction aims mainly to describe how each method behaves, with system-level impacts 
discussed in a later section. 

Shown notionally in Figure 68, surface cooling by means of cooled gas through a porous 
wall is known by different terms in the literature. Distributed film cooling is where a 
normally solid wall is drilled with an array of holes, inline or staggered, using mechanical, 
laser, or water drilling. Transpiration (and more uncommonly, effusive) cooling is where 
the wall is composed of a material that is manufactured to be porous. Examples of such 
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materials include carbon/carbon composites (C/C) [80, 81], ceramic matrix composites 
(CMC) [82], and porous stainless steel manufactured by layering and sintering stainless 
steel mesh [77]. For the purposes of this study, no practical distinction is made between 
distributed film cooling and transpiration — both will be referred to under the umbrella term 
“transpiration.” A transpiration cooling concept cools the surface in two ways. The first 
process occurs as the cool gas travels from a subsurface plenum to the surface through pores 
in the wall. As the gas travels through the wall, the material is cooled convectively as heat 
travels from the wall to the fluid. The second process occurs once the fluid emerges from 
the surface. If the mass flux ratio between the injectant and the main flow is high enough, 
a protective, cool layer of gas will form on the wall and reduce the heat load traveling from 
the previously hot shear layer to the wall. This process is nominally called “film cooling”, 
not to be confused with the film cooling method of introducing coolant to the flow. 



Film cooling as a primary cooling method is typically achieved by injecting cool fluid 
tangentially or at an angle to the wall. Both of these methods are shown notionally in Figure 
69. This flow may be introduced via a slot or ports on a backward- facing step, although 
transverse injection is another possible, albeit less conventional, method of injection. The 
effectiveness of him cooling is rooted in the “him cooling” process mentioned previously: 
a layer of gas provides wall protection from a high heat load that would otherwise be 
introduced by a hot, upstream shear layer. This protection is introduced not only by having 
a cool layer in the vicinity of the wall, but also by altering the structure and thickening 
the thermal boundary layer as the cool gas mixes with the upstream shear layer, thereby 
reducing heat loads and temperatures present in the boundary layer. 

Most of the literature, reflected in the general description of the cooling methods to this 
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Figure 69: Notional diagram of film cooling via a slot for both (a) tangential and (b) angled 
injection. 

point, has a single-objective focus on the desire to lower heat loads at and traveling through 
the wall. In essence, temperatures experienced away from the wall (as long as they do not 
affect the wall heating) are of little consequence. For example, a film or transpiration cooler 
that uses a combustible gas such as hydrogen as its coolant performs well if the hydrogen 
is successfully able to keep the surface cool even if the edge of the cooling layer begins to 
combust due to the high temperatures in that region [83]. Similarly, if the primary flow 
is composed with a fuel, burning on the coolant edge is also acceptable as long as heat 
loads remain low at the wall. This is because for a cooling concept occurring internally in 
a combustion chamber additional burning occurring away from the wall is not necessarily a 
detriment (and possibly is a benefit) — a combustor is not a region where ignition mitigation 
is desired. For a forebody flow, however, this is a disastrous outcome because of the increased 
wall pressure and associated forebody drag. 

This chapter focuses on the study of autoignition mitigation in the locality of the fore- 
body fuel injectors studied in Chapter 2. Reacting simulations of the mixing regions of 


80 


the stinger, ramp, and strut injectors were evaluated in order to determine the autoigni- 
tion characteristics of each injector. After identifying injectors that showed promise for 
autoignition mitigation, surface cooling was added to the simulations in order to determine 
the parametric effect of coolant species and secondary injection parameters on the presence 
of autoignition in the wall region. 

3.1 Model Selection and Validation 

3.1.1 Chemical Model Selection via Ignition Delay 

The prediction of chemical induction and subsequent ignition is highly sensitive to the chem- 
ical model used in the numerical analysis. Because the ignition delay time has significant 
consequence concerning forebody autoignition, it is critical that the correct chemical model 
be used for the regime of interest. 

The chemistry solver and model used in this study were chosen using the following 
methodology: 

1. Compare the ignition time predicted between solvers. Choose the best solver depend- 
ing on its solution time and its performance against validation data. 

2. Using a selected solver, compare the ignition time predicted between various models. 
Pick the model that best matches validation data. 

3. Using a selected solver and full chemical model, compare the ignition time predicted 
using reduced models. 

There are three solvers whose performance and predictive capability were be compared: 
Cantera, CVODE, and a — QSS. The selection of the chemical solver has importance for 
the construction of the QUIVR code used in Chapter 5 as VULCAN has its own chemistry 
solver, validated in Section 3.1.2. The evaluated solvers were: 

• Cantera 1 , an open source, Python-based chemical kinetics package led in develop- 
ment by Dr. D. Goodwin of Caltech. The package is designed to interface with C++, 
Fortran, Python, or MATLAB. 

1 Code and documentation are obtainable from http://code.google.eom/p/cantera/ 
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• The CVODE package of SUNDIALS 2 , an open source nonlinear and differential 
equation solver suite developed by Lawrence Livermore National Laboratories. The 
package is able to interface with C++, Fortran, or MATLAB. 

• The a — QSS algorithm, developed by Mott et al. [84], which is specifically tuned for 
the solution of stiff chemical equations. 

The Cantera solver was called from a lightweight Python script coded by the author. The 
CVODE and a — QSS routines were called from a Fortran 95 zero-D reactor program also 
coded by the author which can accept any general set species and reactions without pressure 
dependent terms. When comparing chemical models against each other and experimental 
data, the combustion problem was modeled as occurring in a constant pressure reactor. 

The literature indicates at least four ways to determine when combustion occurs, and 
therefore ignition delay time. The first is by locating the time of maximum rate of tem- 
perature change. The second is by finding the maximum rate of OH production [85]. The 
third is by defining a threshold of OH concentration, for example 0.005 mol /m 3 [86]. A 
fourth method defines a threshold ignition temperature [86] . The third definition was used 
in this study because it is also a convenient means to determine ignition location in a 
multidimensional flow simulation. 

For the purposes of comparing solvers, the full Jachimowski chemical model was chosen 
as the baseline chemical model. The Jachimowski model has been used extensively in past 
numerical studies and has been validated over a wide range of conditions [87]. Ignition 
delay times using the model were simulated using a constant-volume combustion code by 
McLain and Rao [88] and compared against the experimental results of Slack and Grillo [85] 
with most simulated delay times being within 20% of experiment. Burning velocities were 
compared against data compiled by Warnatz [89] and Milton and Keck [90] using a premixed, 
ID flame code written by Kendall and Kelly [91]. Generally good agreement between 
simulation and experiment were found both in terms of maximum burning velocity as well as 
burning velocities with varying equivalence ratio. Finally, scramjet combustor simulations 

2 Code and documentation are obtainable from https://computation.llnl.gov/casc/sundials/main.html 
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were conducted at flight enthalpies ranging from Mach 8 to 25 which demonstrated the 
importance of including HO 2 chemistry in the chemical model due to the high sensitivity 
of the results to that chemical species. It was also observed that at flight enthalpies above 
Mach 12, nitrogen chemistry became important. 

The performance of each solver over a fixed temperature range at a pressure of 2.0 
atm is shown in Figure 70. The shock-tube experiments of Slack and Grillo [85] at the 
same ambient conditions were plotted in order to provide a means of comparison. At high 
pressures these experiments provide high quality results, although the quality diminishes at 
lower pressures due to the influence of the boundary layer in the shock tube. The general 
trend is observed to have long ignition delay times at low temperatures with the delay time 
rapidly shortening as temperature was increased. As temperatures increased past the knee 
in the curve, the decrease in ignition delay time becomes more gradual. It is evident that 
regardless of the solver used, the ignition time predictive capability was the same. Because 
there is computational overhead required in using the Sundials and Cantera packages, versus 
a — QSS which can be simply coded, a — QSS was chosen to be the best algorithm for this 
study for the QUIVR code. 

With the solver comparison complete, the effect of using different models was ascer- 
tained. A number of models were chosen, including the Jachimowski, O Conaire [92], 
Rogers and Schexnayder [86], and San Diego [93] models. Appendix D gives a summary 
of each mechanism. Figure 71 shows each model’s predicted ignition times compared to 
the experiment of Slack and Grillo at a pressure of 2 atm. It is observed that each model 
gave a reasonable representation of experimentally determined ignition delay times. The 
largest differences occurred at low initial temperatures. Figure 72 compares the same mod- 
els against an experiment conducted by Bhaskaran [94] at a pressure of 2.5 atm. Generally 
good agreement is again achieved, although the San Diego model agreed closest with ex- 
perimental data at low temperatures. 

Because of the simplicity of the model, its established use in flows with scramjet-like 
conditions, and generally good agreement with experiment, the Jachimowski model was 
chosen to be the model used in the present chemically reacting flow studies. Because 
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Figure 70: Comparing ignition delay time as a function of initial reactor temperature for a 
constant pressure of 2.0 atm using different solvers. 


p = 2.0 atm 



Figure 71: Comparing ignition delay time as a function of initial reactor temperature for a 
constant pressure of 2.0 atm using different chemical models. 
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Figure 72: Comparing ignition delay time as a function of initial reactor temperature for a 
constant pressure of 2.5 atm using different chemical models. 

running the full chemistry model had a significant adverse impact on the computational 
resources required to converge a numerical simulation, a study of the effect of reducing 
the model is undertaken here. This was accomplished by removing selected species from 
the mechanism and related reactions and evaluating the effect on ignition delay time. The 
result of this is shown in Figure 73, which shows the ignition delay time predicted by 
reduced Jachimowski models compared to the 2 atm shock tube experiment of Slack and 
Grillo. Here it is observed that removing all of the nitrogen chemistry as well as hydrogen 
peroxide (H 2 O 2 ) chemistry had a negligible effect on the prediction of ignition delay time 
at these conditions. The relative unimportance of hydrogen peroxide for high temperature 
combustion has been established in the literature [95]. Therefore, the 8-species, 13-reaction 
reduced Jachimowski chemical mechanism was chosen to be the preferred chemical kinetic 
model moving forward. 
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Figure 73: Comparing ignition delay time as a function of initial reactor temperature 
for a constant pressure of 2.0 atm using different reduced models from the Jachimowski 
mechanism. 

3.1.2 Blunt-Body, Premixed, Shock-Induced Combustion 

The VULCAN CFD code in a chemically reacting mode of operation can be validated 
against blunt-body, premixed, shock-induced combustion experiments conducted by Ruegg 
and Dorsey, and H.F. Lehr. In both cases, a spherically-blunted projectile was fired into 
quiescent hydrogen-air at various combinations of Mach number and ambient pressure. 
Computational validation is assumed to be achieved if simulations are able to produce the 
physical mechanisms observed in experiment as well as quantitative matching of different 
experimental quantities such as standoff distances and reaction front frequency (for the 
unstable cases). Therefore the validation questions here are summarized as: 

Can VULCAN... 

• ...successfully reproduce the stability conditions observed in experiment at the same 
conditions? 
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• ...accurately predict the reaction and shock front standoff distances observed in ex- 
periments at the same conditions? 

• ...accurately reproduce the same reaction front frequencies as observed in experiments 
at the same conditions? 

These validation questions were chosen due to the importance of ignition delay time and 
chemical model on the behavior of autoignition. 

The computational grid used to model the spherically-blunted experimental projectiles 
is shown in Figure 74. Although the important physical mechanisms responsible for the 
unsteady behavior observed in experiment was hypothesized to occur in the stagnation 
region, the flow is simulated up until the shoulder region of the projectile. The grid was 
non-dimensionalized to the sphere radius in order to allow the grid to be sized easily to 
different projectile radii. This also produced the final results with dimensions normalized 
to sphere radius, a typical normalization used in such flows. For each case the grid had 
800 points normal to and 256 points along the projectile surface. The fineness of the grid 
normal to the surface allowed the resolution of the mechanisms responsible for unstable 
behavior in the reaction zone. All cases were simulated with the Jachimowski 8-species, 13- 
reaction model [87] as this model was concluded to be sufficient for the accurate prediction 
of ignition delay times. During the chemically reacting integration the global CFL number 
was set to be 5.0, which corresponded roughly with a time step of 7 ns. 

In their ballistic experiments, Ruegg and Dorsey [96] shot 20 mm diameter, spherical 
projectiles into quiescent hydrogen-air and measured, among other things, the shock- and 
reaction-front standoff distances at different Mach numbers and ambient pressures. They 
also published Schlieren images that clearly showed different stability regimes for various 
combinations of Mach number and ambient pressure. They linked the relative standoff 
distance between the shock and reaction front to ignition delay time and made estimates 
of that time quantity for each experimental case. Both qualitative and quantitative com- 
parisons can be made between VULCAN and experiment. First, qualitative comparison 
between the simulated reacting flowfield and experimental shadowgraphs can be made to 
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Figure 74: The grid used in the validation study is shown with every tenth grid point 
normal to the surface and every third grid point tangential to surface. 

ensure the same same stability regime is being reproduced at the same conditions. Sec- 
ond, quantitative comparison can be made between the computed reaction front and shock 
standoff distances and those determined from experimental Schlieren imagery. 

Qualitative comparisons of VULCAN results were made against experimental imagery. 
Conditions were chosen from Table 12 in order to reproduce the three major regimes of 
stability observed in experiment. Figure 75 shows the VULCAN density solution of case 
RD-1 compared to the Schlieren image captured at the same experimental conditions. This 
case demonstrates steady behavior of both the reaction front and the bow shock. Such a 
case is characterized by smooth reaction and shock fronts with no direct interaction between 
either except for the temperature being raised behind the shock such that chemical induction 
may occur. Qualitatively the numerical results matched the experimental results as the 
Schlieren image also indicated stable reaction and shock fronts. 

The regular, unsteady behavior characteristic of case RD-2 is observed in Figure 76. 
Here the shock and reaction fronts exhibited periodic motion normal to the surface of the 
nose of the projectile with the distortions then being convected downstream. The greatest 
disturbance was observed in the reaction front. The density flowfield shows all of the salient 
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Figure 75: (a) Numerical shadowgraph of a Mach 5, 20 mm spherically blunted projec- 
tile fired into 0.1 atm quiescent H2-Air is shown and compared to the (b) experimental 
shadowgraph [96]. 
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Table 12: Cases simulated for the Ruegg and Dorsey experiments are listed. All conditions 
are at 300 K with a spherical diameter of 20 mm 


Case ID 

Mach 

p (atm) 

Stability Regime 

RD-1 

5.0 

0.1 

Steady 

RD-2 

4.9 

0.25 

Regular, Unsteady 

RD-3 

4.8 

0.5 

Large Disturbance, Unsteady 


features which encompass the mechanism of stability including contact discontinuities con- 
necting each corrugation to the bow shock. This mechanism is discussed in more detail in 
Chapter 5. 

The final regime of stability, the large disturbance regime, is reproduced in case RD- 
3 as shown in Figure 77. The experimental results at the same conditions show that a 
large amplitude disturbance was present on the reaction front. Numerical results also show 
this behavior, which was due to a deflagration-to-detonation transition (DDT) that caused 
the reaction front to advance upstream in a superdetonation and interact directly with the 
shock front. It is difficult to discern other details from the Schlieren image, but qualitatively 
both flowfields have the same large amplitude behavior with respect to the reaction front. 
Therefore there is qualitative agreement between the numerical simulation and experiment, 
lending confidence that VULCAN is able to simulate the salient physics of this regime, in 
this case the presence of DDT. 

Attention is now given to a quantitative comparison of VULCAN’s predictive capability 
with results from the experiments. For each combination of projectile Mach number and 
ambient pressure, Ruegg and Dorsey noted the shock and reacting front standoff distances. 
The data are shown in Figures 78 through 80 and present the experimental and simulated 
standoff distances (normalized with respect to projectile radius) for the three cases presented 
up to this point. Experimental trendlines are also indicated. Some experimental values lying 
far from the trendline are actually from a separate experiment from the Ruegg and Dorsey 
experiments that used a different method of visualization. Therefore, the trendlines ignore 
the effect of these outlier cases. The extracted standoff distances from the simulated cases 
lie comfortably within the experimental scatter and near the fitted trendline. This shows 
that VULCAN is accurately predicting the quantitative shock and reaction front standoff 


90 



(b) 


Figure 76: (a) Numerical shadowgraph of a Mach 4.9, 20 mm spherically blunted projec- 
tile fired into 0.25 atm quiescent H2-Air is shown and compared to the (a) experimental 
shadowgraph [96]. 


91 



Bow Shock 



Reaction 
Front — 


; 

lips 


■ 

(b) 


Figure 77: (a) Numerical shadowgraph solution of a Mach 4.8, 20 mm spherically blunted 
projectile fired into 0.5 atm quiescent H2-Air is shown and compared to the (b) experimental 
shadowgraph [96]. 
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Figure 78: A comparison is shown between shock and reaction front standoff distances for 
experiment and the present simulation at 0.1 atm. 


distances. 

The final validation question relates to the ability of VULCAN to predict the reaction 
front frequency for a regular, unsteady regime case. If the reaction front is positioned far 
enough away from the projectile nose, the periodic reaction front behavior can be considered 
to be inviscid in nature, with the induction time between the shock and reaction front being 
the dominant parameter of interest. Hartmuth F. Lehr conducted early experiments [10] of 
reaction front frequency ahead of 15 mm diameter projectiles fired into 0.421 atm, 292 K, 
quiescent hydrogen- air. 

The frequency results presented herein used the same computational grid and models 
used in the Ruegg and Dorsey simulations. The only differences are the ambient conditions 
and projectile diameter. The solution was considered to be converged whenever the tran- 
sients passed from the system, leaving only the physical, periodic behaviors. The converged 
Mach 4.47 flowfield in Figure 81 shows a numerical shadowgraph with both the reaction 
and shock fronts clearly visible. The reaction front shows a corrugated shape indicative of 
the regular, unsteady regime. 
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Figure 79: A comparison is shown between shock and reaction front standoff distances for 
experiment and the present simulation at 0.25 atm. 



Figure 80: A comparison is shown between shock and reaction front standoff distances for 
experiment and the present simulation at 0.5 atm. 
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Figure 81: Numerical shadowgraph for the Mach 4.47 projectile case with stoichiometric 
H2-air at 0.421 atm and 292 K ambient conditions. 

Further interrogation of the flowheld was required in order to obtain the reaction front 
frequency. Lehr calculated his reaction front frequencies by simply measuring the spaces 
between corrugations in the wake of the projectile and combining that with an estimate of 
post-shock velocity. The current methodology simply interrogated the stagnation streamline 
at defined time intervals, a density plot of which is shown in Figure 82. It is observed that 
the regular regime behavior closely matches the model shown in Figure 122. 

Once the stagnation streamline distance versus time plots were obtained, it was a rela- 
tively straightforward process to obtain reaction front frequencies. The process was to take 
a time- varying slice at a chosen value of n and then compute a discrete Fourier transform to 
obtain the frequencies present at this slice. Any arbitrary location can be chosen behind the 
reaction front, and in this case it was chosen to take the slice at the wall. The time- varying 
water mass fraction at the wall for the Mach 4.47 case is shown in Figure 83. Note that 
in addition to the high frequency, high amplitude behavior there was also a low frequency, 
low amplitude behavior superimposed on the trend. This was a physical behavior as this 
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Figure 82: Density contours on the stagnation streamline for the Mach 4.47 case are pre- 
sented as a function of distance from the projectile nose and time. 

case resides in what is referred to as the Envelope Regime [97], where two frequencies are 
present in the flow. The low frequency behavior is a gradual advance and retreat of the 
reaction front as a function of time and would not be visible on a Schlieren image as it does 
not manifest itself by a density change as with the high frequency behavior. 

A discrete Fourier transform of the Mach 4.47 data shown in Figure 83 with a rectangular 
windowing function is plotted in Figure 84. There is a low-frequency peak at 87 kHz and a 
high frequency peak at 434 kHz. Note that there is some variance to these peaks because of 
uncertainty due to spectral leakage. At the same condition, Lehr computed a reaction front 
frequency of 425 kHz, resulting in a 2.1% relative error between the numerical simulation 
and experiment. A lower frequency is not reported by Lehr, but numerical work by Myles 
Sussman [98] reports a value of 74 kHz, resulting in an 8% relative error with respect to 
the mean. 

It has been shown that VULCAN is able to accurately predict the characteristics of 
time varying, chemically reacting flow fields. The selection of premixed, shock-induced 
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Figure 83: Water mass fraction as a function of time at the projectile wall for the Mach 
4.47 case. 



Figure 84: The discrete Fourier transform taken of the wall water mass fraction for the 
Mach 4.47 case. 
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combustion as a method of validation was a good one because the problem contains many 
different mechanisms against which to evaluate a code. First, depending on the regime 
studied, the reaction front varies with time with a prescribed frequency. The ability of the 
code to capture this behavior gives credibility in its ability to properly simulate the physical 
dynamics occurring in the flow field. Second, the reaction front location and frequency are 
strongly dependent on the ignition time of the fuel-air mixture behind the shock wave. 
Neither can be accurately determined without a proper handling of the chemistry and 
careful selection of the chemical model. Because the VULCAN code was able to predict the 
regimes of stability, the shock and reaction front standoff distances, and the reaction front 
frequencies as compared to experiment, it is concluded that VULCAN and the Jachimowski 
reduced chemical model are valid for the regimes of interest studied in this dissertation with 
the assumption that the methods are valid and extrapolate to higher Mach numbers. 

3.2 Flow Simulation and Metrics of Performance 

3.2.1 Flow Conditions and Computational Domain 

Unless otherwise indicated, the flow conditions were the same as those used in the mixing 
study of Chapter 2. For the initial autoignition studies, the same grid topology was used as 
those used in the mixing study. Transpiration cooling cases used the same grid topology as 
the uncooled geometry with the change being that transpiration was applied as a boundary 
condition. This also holds true for angled film cooling. For axial film cooling, a backward- 
facing step was placed at the wall at the cooling location with the wall-normal grid spacing 
held constant and equal to the grid spacing at the wall without the backward facing step. 

3.2.2 Measurement of Oxygen Burning 

In this study combustion efficiency and the change of oxygen burned in the sampling region 
were the means of quantitatively gauging the effectiveness of any given autoignition mitiga- 
tion strategy. Combustion efficiency is a numerical value between zero (no combustion) and 
one (complete combustion where all reactants not in excess are consumed) that is typically 
employed in a fashion where higher values occurring sooner are better. Because the current 
study was situated on the vehicle forebody, however, it was desired to delay combustion as 


98 



Figure 85: Top -level schematic of performance sampling region that spills oxygen (for a 
system where fuel is in excess). 

long as possible, thereby making combustion efficiency a metric to be minimized. There are 
many ways to define combustion efficiency depending on the particular application. Typical 
measures of combustion efficiency include measures of 


1. Water produced in the domain of interest compared to the reactant not in excess. 

2. Consumed reactant mass not in excess compared to the initial reactant mass. 

3. Current temperature compared to the adiabatic flame temperature. 

In high speed, airbreathing propulsion systems combustion typically occurs in an en- 
closed space (e.g. combustion chamber and/or nozzle). Because of this, combustion effi- 
ciency metrics generally do not account for spillage of reactants such as that shown notion- 
ally in Figure 85. Therefore, a new combustion efficiency metric that takes reactant spillage 
into account will now be proposed. 

For the performance sampling region notionally shown in Figure 85, reactants entering 
the system are spilled in some proportion relative to the inflow. Because hydrogen is 
typically in excess for the flows studied here, oxygen was the tracked reactant. Figure 86 
shows a notional plot of oxygen mass flow as a function of distance along the body for both 
a reacting and nonreacting simulation. The two main quantities of interest for computing 
combustion efficiency, nonreacting and burned mass flow rates, are indicated in the figure. 
Note that in order to obtain the burned oxygen mass flow rate, both a nonreacting and 
reacting simulation of the same flow are required. 

Once the nonreacting and reacting mass flow rates are obtained for a particular case, it 
is relatively straightforward to define a combustion efficiency parameter that satisfies the 
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Figure 86: Notional oxygen mass flow rates as a function of distance along sampling region 
with spilled and burned quantities indicated. 

requirements that it be equal to zero when no burning has occurred at a particular location 
and equal to unity when all available oxygen at a particular flow station has been consumed. 
The combustion efficiency metric that takes into account reactant spillage is defined as 

Vc = = 1 - ( 8 ) 

17102, nonreac ,nonreac 

To obtain the combustion efficiency at a particular station, the reactant mass flow rate for 
the nonreacting and reacting simulation must be obtained. For codes that do not output 
these quantities automatically, this process is the summation of cellular mass flow rates 
multiplied by the mass fraction of the reactant of interest in that cell. 

For comparing different geometries, a more general metric was used, the change in 
oxygen burned in the sampling region. Different cooling methods and secondary injec- 
tion conditions will change the spillage properties of a particular geometry being analyzed. 
Therefore, spillage must be taken into effect directly in the metric measuring oxygen burned. 
This equation, using Figure 85 as a reference, is 

^ 02 , burned ( r i rl 02, spill + hlQjjCmt) (9) 

The change in oxygen burned can then be computed by evaluating this expression sepa- 
rately for the baseline reacting geometry and for the cooled geometry being investigated, 
subtracting the former from the latter. 
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3.3 Baseline Simulations of Autoignition 


For the present results, the geometries shown in Figures 37, 38, and 39 were simulated at 
the conditions in Table 3. The simulations were allowed to first converge without chemistry, 
at which point chemistry was activated and the simulation was allowed to reach convergence 
again. Upon convergence, chemical ignition is tracked by the water production rate in units 
of kg/m 3 /s. Because the fuel is pure hydrogen, the production of water was also tracked in 
the flow. 

The converged solution of the stinger flush- wall injector showed the greatest amount 
of chemical reaction of the three injector concepts tested. Figure 87 shows centerline and 
cross-stream slices of water production and water mass fraction. The water production are 
good indications of where ignition begins. Two primary ignition points are evident: behind 
the bow shock of the injector plume and laterally on the centerline between neighboring 
injectors. The ignition region behind the bow shock appears to contribute the most to the 
combustion occurring in the flow by wrapping the entire fuel plume in shock-heated air 
and easily causing ignition as hydrogen diffuses out of the fuel rich plume. This shock- 
layer burning has been observed in OH-PLIF experimentation by Ben-Yakar, et al. [99] 
for circular jet injection into a hypervelocity, Mach 13 flight enthalpy flow. The stinger 
injector geometry, however, lacks the upstream recirculation burning observed with circular 
injection due to the sharp leading edge of the injector eliminating such a recirculation zone. 
It is unclear how a shock layer ignition zone can be mitigated using a wall cooling strategy 
because this ignition zone is caused by high temperatures behind the bow shock, not high 
temperatures in the boundary layer. 

The secondary region of ignition on the symmetry plane between neighboring injectors 
was due to high vorticity that entrained fuel and mixed it with the hot boundary layer near 
the symmetry plane. This created a large, high-temperature region in region of the flow 
that also was a secondary location for autoignition. Such an autoignition region may be 
mitigated by wall transpiration in the locality of the symmetry plane. 

A plot of combustion efficiency as a function of distance bethind the stinger injector 
is shown in Figure 88. In this plot the amount of combustion throughout the flowfield 
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Figure 87: Stinger injector contours of (a) water production and (b) water mass fraction 
in the flow domain of interest. Flow is from the upper left to the lower right. 
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Figure 88: Mass flow rates of oxygen and combustion efficiency downstream of stinger 
without autoignition mitigation. 

was quantified and accounts for spillage in the computation. It is observed that by the 
time the outflow plane was reached, the combustion efficiency achieved a value of 45%, 
indicating extensive combustion and fuel depletion. While this may be a good trend for 
internal combustion, it is bad for forebody injection. Because shock-heating of the air plays 
an important part in the trend, him cooling may be an infeasible approach for autoignition 
mitigation. The high losses experienced by the transverse injector in a previous study [100] 
combined with the present autoignition characteristics disqualifies the transverse injector 
as a feasible method for forebody injection. 

The strut injector experienced a moderate amount of autoignition in the wake of the 
concept. Two distinct zones of ignition are observed, one away from the wall in the fuel 
plume vortex and another at the wall, both shown in Figure 89. Shock interactions behind 
the strut cause subsequent ignition which is extinguished relatively rapidly as the products 
are folded into the cooler fuel plume. Downstream of the first ignition point, a second 
ignition occurred at the wall as the fuel plume entered into the hot boundary layer region. 
Downstream of the second ignition point combustion is sustained near the wall, evident 
by the continuous production of OH at the wall. While the ignition region at the wall 
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may be addressed by wall cooling, the small ignition region in the fuel plume may require 
modification to the injector geometry in order for complete autoignition mitigation to occur. 
The extent of the combustion occurring for the strut concept is shown quantitatively in 
Figure 90. Much lower levels of combustion are evident versus the stinger concept with 
combustion efficiency approaching only around 12%. 

The unswept ramp injector with wedge feature experiences a single region of autoignition 
downstream of the injector at the wall and near the outflow plane. Contours of water 
production rate and water mass fraction are shown in Figure 91. Here it is again seen 
that ignition occurred near a region where the fuel entered the hot boundary layer near the 
flow domain exit. Because this is the only zone of autoignition in the flow and it is near 
the wall, wall cooling is a strong candidate for mitigating ignition in this area. Figure 92 
quantitatively confirms that no substantial combustion is occurring in the flow, with the 
combustion efficiency trend mostly reflecting numerical errors. 

Through investigation of autoignition for the stinger, strut, and unswept ramp injector 
with base wedge, the strut injector distinguishes itself as an attractive candidate for fur- 
ther study in autoignition mitigation. The stinger injector is not an attractive candidate 
for further study simply due to the large extent of combustion present, its difficulty for 
mitigation, and other losses incurred even without reaction present. Similarly, the ramp 
injector is not an attractive candidate to study due to its general redundancy with the 
strut case. It, like the strut, has autoignition occur once the fuel plume enters into the hot 
boundary layer. Its main distinction is that this effect occurs farther downstream than the 
strut injector. It can be assumed that strategies that are successful for the strut injector 
will also be successful for the ramp injector. Therefore, the study of autoignition mitigation 
strategies will concentrate on the strut injector from this point onward. 

3-4 Autoignition Mitigation Strategies for Forebody Flows 
3.4.1 Overview of Wall Cooling Strategies 

Wall cooling using both wall transpiration and film cooling were evaluated in this study to 
ascertain their performance in mitigating forebody autoignition. In both cases a secondary 
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Figure 89: Strut injector contours of (a) water production and (b) water mass fraction in 
the flow domain of interest. Flow is from the upper left to the lower right. 
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Figure 90: Mass flow rates of oxygen and combustion efficiency downstream of strut without 
autoignition mitigation. 

injectant was introduced into the main flow at the wall in order to lower the temper- 
atures encountered in the boundary layer below those that would cause autoignition to 
occur. The physical modeling in this study makes no distinction between the transpiration 
and distributed film cooling strategies, with the main distinctions arising due to material, 
structure, and subsystem design. Both strategies will henceforth be referred collectively as 
“transpiration cooling.” Transpiration is shown notionally in Figure 68 and film cooling is 
shown notionally in Figure 69. 

Both transpiration and film cooling have benefits and detriments beyond the system- 
level metrics shown in the table. For example, transpiration is typically a good strategy for 
providing localized cooling at the wall throughout the region transpiring coolant. Beyond 
the transpiration patch, however, cooling capacity drops off quickly. Wall transpiration 
also typically provides good cooling effectiveness while also minimizing disturbance to the 
main flow [28, 78]. Film cooling, which injects fuel at a single location in order to cool the 
wall downstream, typically requires more mass flow from the injector to achieve the same 
amount of wall cooling due to decreased cooling capacity as the injectant flows downstream. 

A blowing parameter is defined in order to control the amount of secondary gas flowing 
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Figure 91: Unswept ramp injector contours of (a) water production and (b) water mass 
fraction in the flow domain of interest. Flow is from the upper left to the lower right. 
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Figure 92: Mass flow rates of oxygen and combustion efficiency downstream of strut without 
autoignition mitigation. 

into the wall boundary layer. This parameter is a function of the mass flow rate both 
through the primary injector and through the secondary injector and is defined as 

B = mc °° l . (10) 

dicooi T TTlinj 

where m coo i is the mass flow rate issuing from the film or transpiration cooler and rhi n j is 
the mass flow rate issuing from the injector. 

The secondary gases in this study were chosen to be hydrogen, helium, and nitrogen. 
Hydrogen was considered because it is the same gas as the fuel and, therefore, would only 
require the system to provide a single fuel storage solution. It also has a high heat capacity 
by mass, which will possibly provide for lower temperatures as the transpired gas mixes 
with the main flow. A detriment to using hydrogen is its chemical activity as it readily and 
rapidly burns. 

Helium and nitrogen provide respectively low and high molecular weight alternatives to 
hydrogen as secondary gases. Their primary benefit compared to the use of hydrogen is 
that they are both inert at the temperatures of interest in the current study. Both have 
detriments compared to hydrogen in that they would require separate storage and cryogenic 
cooling requirements from the fueling system. In addition, nitrogen has poor heat capacity 
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by mass compared to hydrogen, although it is a inexpensive and prevalent. Helium can 
provide a comparable heat capacity by mass to hydrogen, although it is more expensive and 
less available than nitrogen. 

Other examples of coolant candidates include water, ammonia (NH 3 ), and methane 
(CH 4 ). Water was initially excluded due to risk of oxidation of cooling passages, ammonia 
due to the increased cooling mass required, and methane due to the possible presence of 
graphite coking in the potentially hot cooling passages [ 101 ]. 

The selection of the blowing parameter value has a direct impact on the global equiva- 
lence ratio of the system, with higher blowing parameters necessarily increasing the global 
equivalence ratio and mass flow rate of hydrogen demanded on the system. Because sec- 
ondary gas is intended to inhibit autoignition, regardless of whether the selected gas may 
burn or not, any increase in fuel mass flow rate can be construed to be parasitic on the 
overall system performance. Figure 93 shows the effect on the global equivalence ratio by 
increasing the blowing parameter as given by the equation 

rilinj OF 

G rhair (1 - B)Yox 

where OF is the oxidizer-to-fuel mass ratio and Vox is the stoichiometric mass fraction 
of the oxidizer. This shows that, for a general fuel-air system, as blowing parameter is 
increased, the total global equivalence ratio supplied by the injection and cooling systems 
increases. 

Because the mass flow rate of fuel — coolant gas also being counted as fuel regardless 
of flammability — factors directly into the calculation of specific impulse, a system-level 
efficiency metric, the effect of blowing parameter on the change in this parameter can be 
assessed. Figure 94 shows the relative change in specific impulse on the choice of blowing 
parameter, assuming the thrust produced is unaffected. Here it is observed that low blowing 
parameters can have a rather drastic effect on the performance of a system using wall 
cooling. Therefore, the importance of using only the minimum amount of secondary gas is 
emphasized. 
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Figure 93: Global equivalence ratio as a function of blowing parameter 



Figure 94: Relative change in specific impulse as a function of blowing parameter 
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3.4.2 Wall Transpiration 


Wall transpiration (i.e. blowing) of cold gas into the hydrogen/air main stream was chosen 
as one method for autoignition mitigation in this study. Transpiration has been used in the 
past for wall cooling in Mach 2-3 flows with flush-wall injection. While it is not the only 
way to provide cooling near the wall, transpiration has good wall cooling effectiveness while 
also minimizing the disturbance to the main flow. 

For the sake of computational simplicity, fluid flow through transpiration injection ports 
was modeled using a representative bleed model rather than gridding and simulating distinct 
wall holes and a plenum. This allowed rapid evaluation of the effect of different porosities 
and mass flows without altering the grid topology. Note that the use of such a model is 
particularly attractive for blowing studies because of the insensitivity of the solution to the 
shear stress at the wall for this case [102]. The ability of VULCAN to accurately model 
this phenomenon has been studied by Baurle and Norris [103]. 

The transpiration conditions used in the study are given in Table 13 for the transpiration 
schematic shown in Figure 95. The transpiration model parameters in VULCAN are the 
mass flow through the surface, the total temperature of the fluid in the plenum, and the 
porosity of the wall. The porosity is defined as the ratio of the combined area of the 
transpiration holes in the wall to the total area of the transpiration wall. The porosity 
was estimated from past experiments by Byington et al. [28]. Therefore, the porosity for 
this study was set to be 2% in order to create the greatest possible plenum pressure. It 
is also the lower limit of validation for the transpiration model. The total temperature of 
the transpirant was set to be near the wall temperature. Besides the species and blowing 
parameter, the size of the transpiration patch was another variable for the current study. 


Table 13: Transpiration array conditions 


Blowing parameter ( B ) 

5.0 - 20.0% 

Wall porosity (P) 

2.0% 

Transpirant total temperature 

500 K 

Wall temperature 

500 K 

Species 

H 2 , He, and N 2 
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Figure 95: Schematic diagram of wall transpiration. 


3.4.3 Film Cooling 

In addition to transpiration, film cooling was another candidate method for autoignition 
mitigation. Film cooling is distinct from transpiration in that the injectant is supplied 
via a channel or row of injectors instead of through an array of perforations in the wall. 
Furthermore, if the flow is to be supersonic, a contoured nozzle would be present in the film 
cooling channels. While injection through a channel transverse to the main flow is possible, 
axial injection is preferable due to the desire to keep the fluid near the wall. The use of 
tangential injection also reaps minor thrust benefits to the vehicle. Furthermore, when a 
low molecular weight gas like hydrogen is injected at the wall, reductions in skin friction 
may be achieved [104], 

The mass flow rate supplied to the film cooler, like the transpiration case, was parametrized 
in terms of the blowing parameter defined in Equation 10. Other parameters to be studied 
include the species injected and the exit temperature of the film cooling channel, summa- 
rized in Table 14 for the schematics shown in Figure 96. 


Table 14: Film cooling conditions 


Blowing parameter ( B ) 

5.0 - 20.0% 

Channel height 

5 mm 

Injection angle 

0°, 10°, and 45° 

Exit temperature 

100, 740 K 

Species 

H 2 , He, and N 2 


While some film cooling concepts inject gas transverse to the flow, the current study 
was limited to concepts that contain a substantial tangentially-directed component. Three 
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Figure 96: Schematic diagrams of film cooling via a slot for both (a) tangential and (b) 
angled injection. 

configurations will be tested: parallel injection, 10° injection, and 45° injection. Parallel 
injection is chosen to provide the least disturbance possible to the primary flow, while 
the 10° and 45° injection cases will both necessarily create a shock wave and some flow 
deflection. The reason for testing these angled film cooling concepts was to deflect the hot 
shear layer away from the wall and to see if such benefits offset the creation of a shock wave 
and higher temperatures between the film cooling shock and the secondary flow. 

3.5 Simulation of Transpiration and Film Cooling 
3.5.1 Simulation Matrix 

The performance evaluation of film cooling and transpiration was conducted in a systematic 
fashion to allow the best use of computational and research time. The strategies were 
parametrized in terms of species, blowing parameter, and nozzle exit temperature for film 
cooling and species, blowing parameter, and patch width for transpiration cooling. For the 
fuel mass flow rate utilized in the previous chapter, Table 15 shows the required film cooling 
or transpiration mass flow rate for values of blowing parameter ranging between 5 and 20%. 

For evaluating him cooling cases, the systematic methodology for evaluating each con- 
cept consisted of first computing a pathfinding simulation at a particular combination of 
species and nozzle exit temperature with B = 20%. This was a test case that showed 
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Table 15: Coolant mass flow rates for each value of blowing parameter. 


B[%\ 

TTlfc [kg/s] 

20 

0.003275 

15 

0.002312 

10 

0.001456 

5 

0.00069 


Table 16: Simulation matrix for parallel film cooling cases with sonic injection. Checkmarks 


cases that were run and x-marks inc 

icate cases not run. 

Species 


H 

2 


He 

n 2 

Te [K] 

B [%] 

20 

15 

10 

5 

20 

15 

10 

5 

20 

15 

10 

5 

100 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

740 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 


if the particular combination of species and exit temperature were feasible (i.e. substan- 
tially reduces the presence of autoignition on the forebody). If a particular combination 
of species and exit temperature reduced autoignition substantially, the blowing parameter 
was gradually lowered to evaluate the effect of that parameter. Note that the mitigation of 
autoignition is not “all or none” — one can not hope to employ a strategy that will eliminate 
autoignition at all points downstream of the exit of the film cooler. This is because there will 
always be a distance downstream where the primary and coolant streams have equilibrated 
and the boundary layer once again begins to increase in temperature. The best outcome 
is to delay the presence of ignition until the point where it is desired (i.e. the combustor 
entrance) . 

The film cooling cases tested for this study are shown in Tables 16 through 18, which 
cover the parallel, 10 o, and 45o film cooling strategies. For each case the film cooling channel 
began at the x = 200 mm station. The checkmarks on the tables show both the pathfinding 
cases (at B = 20%) and those additionally run cases to evaluate the effect of decreasing 
blowing parameter. 

The transpiration cases tested for this study are shown in Table 19. Widths are percent- 
ages of the full width of the flowpath centered on the centerline and lengths are percentages 
of the length from the beginning of the transpiration patch at x = 200 mm to the notional 
cowl closure location at x = 713 mm. Similar to the film cooling cases, transpiration cases 
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Table 17: Simulation matrix for 10 o film cooling cases with sonic injection. Checkmarks 


indicate cases that were run and x-marks indicate cases not run. 


Species 


h 2 



He 

n 2 

Te [K] 

B [%] 

20 

15 

10 

5 

20 

15 

10 

5 

20 

15 

10 

5 

100 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

740 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 

/ 


Table 18: Simulation matrix for 45o film cooling cases with sonic injection. Checkmarks 


cases that were run and 

x-marks inc 

icate cases not 

run. 



Species 


h 2 


He 

n 2 

Te [K] B [%] 

20 

15 10 

5 

20 

15 10 

5 

20 

15 10 

5 

100 

/ 

X X 

X 

/ 

/ / 

/ 

/ 

X X 

X 

740 

/ 

X X 

X 

/ 

/ / 

/ 

/ 

X X 

X 


at high blowing levels were tested first in order to establish the effectiveness of the strategy 
at those conditions. If negative net oxygen relative to the baseline was consumed at those 
levels, the blowing parameter was reduced to test the limits of the combination of species 
and patch size. 

3.5.2 Results of Film Cooling 

For each strut case run in Tables 16 through 18, the change in burned oxygen mass flow 
rate with respect to the baseline case was computed. Figure 97 shows a centerline numerical 
shadowgraph superimposed with water mass fractions for the baseline case without any film 
cooling. 

The results are shown in Figure 98 for the case of slot injection of hydrogen. Negative 
values of change in oxygen burned indicate cases where less oxygen was burned in the control 


Table 19: Simulation matrix for transpiration cooling cases with a porosity of 2%. Check- 
marks indicate cases that were run and x-marks indicate cases not run. 
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Injector Wake 
Flow 



x = 150 mm 


Figure 97: Baseline, reacting flow without autoignition mitigation is shown as numerical 
shadowgraph overlaid with oxygen water mass fractions. 

volume relative to the baseline case without any cooling. It was shown to be possible to 
inject hydrogen at the wall and have it mitigate autoignition that would otherwise occur 
without the secondary stream itself autoigniting. This was an important result at the system 
level because it supports the feasibility of a single fuel storage and cooling system for the 
purposes of primary and secondary (coolant) injection. While relatively large values of the 
blowing parameter were required for elimination of wall ignition, this could be addressed 
in a flight design by splitting the fueling between the primary injector and the wall slot so 
that they combine to the desired global equivalence ratio. 

In Figure 98(a) it is shown that in general, for tangential hydrogen coolant injection, 
lower coolant exit temperatures were more advantageous to reducing the amount of oxygen 
burned in the sampling region. This may be related to the increased heat sink capacity of 
the lower temperature coolant. In addition, the percent reduction in oxygen mass burn- 
ing saturates around a 70% reduction when high coolant flow rates are used. Such cases 
succeed in eliminating autoignition at the wall but do not greatly affect the mid-stream 
autoignition zone, as seen in Figure 99. At higher coolant temperatures and blowing rates, 
the secondary hydrogen coolant in Figure 100 exhibited some initial burning while still 
mitigating a large percentage of the wall burning that would otherwise occur. Below 15% 
blowing, the high temperature injectant experienced a rapid drop in performance while the 
drop in performance for the low temperature injectant occured below 10% blowing. These 
drops in performance reflect the presence of wall burning of hydrogen moving upstream as 
the blowing parameter was reduced. At 5% blowing, slightly more oxygen was burned rel- 
ative to the baseline case for the high temperature injectant while for the low temperature 
injectant the amount of oxygen burned nearly broke even. 
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Tangential Hydrogen Injection 
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Te = 740 K 



10 deg Hydrogen Injection 



45 deg Hydrogen Injection 
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Figure 98: Changes in mass flow of oxygen burned for secondary hydrogen injection is 
shown. 
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Injector Wake 
Flow 



Coolant Flow 


x = 150 mm x = 200 mm x = 713 mm 

Figure 99: Numerical shadowgraphs overlaid with water contours Yh 2 o = 0.01, 0.02, ..., 
0.22 for tangential film cooling with hydrogen. Coolant temperature of 100 K and B = 20%. 


Injector Wake 
Flow 



Coolant Flow 


x = 150 mm x = 200 mm x = 713 mm 

Figure 100: Numerical shadowgraphs overlaid with water contours Yh 2 o = 0.01, 0.02, ..., 
0.22 for tangential film cooling with hydrogen. The coolant temperature is 740 K and B = 
20 %. 

Ten degree slot injection performance for the low temperature hydrogen is shown in 
Figure 98(b) and indicates significantly reduced effectiveness of hydrogen coolant. For the 
high temperature coolant at a 20% blowing rate, a 100% increase in oxygen burning was 
realized as combustion occured almost immediately at the location of injection as shown in 
Figure 101(b). This was due to the creation of a shock wave at the leading edge of the film 
cooling slot that caused extensive amounts of coolant burning. Burning was reduced for the 
high temperature coolant by reducing the hydrogen mass flow rate to the slot. This may be 
due to two effects: reduction in strength of the leading-edge shockwave and less combustible 
fuel having been supplied to the wall region. This combustion zone was delayed relative 
to the baseline by lowering the injectant temperature as shown in Figure 101(a), although 
the lowered oxygen burning rates were marginal. These marginal gains by low temperature 
ignition were decreased by lowering the blowing parameter, roughly breaking even at 10% 
blowing, with the amount of burning increasing past 5% blowing. 

Increasing the angle of secondary injection to 45° caused a 180% increase in oxygen 
burned in the injector wake for a high temperature coolant. This can be observed in the 
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(a) 


Injector Wake 
Flow 



Coolant Flow 


x = 150 mm x = 200 mm x = 713 mm 


(b) 


Figure 101: Numerical shadowgraphs overlaid with water contours Yh 2 o = 0.01, 0.02, ..., 
0.22 for 10 degree film cooling with hydrogen. B = 20% and ambient temperatures of (a) 
100 K and (b) 740 K are shown. 


Injector Wake 
Flow 



x = 150 mm x = 200 mm x = 713 mm 


Figure 102: Numerical shadowgraphs overlaid with water contours Yh 2 o = 0.01, 0.02, ..., 
0.22 for 45 degree film cooling with hydrogen. The coolant temperature is 740 K and B = 
20 %. 


bar chart in Figure 98(c) and the numerical shadowgraph in Figure 102. This was primarily 
due to the presence of a stronger shock wave produced by greater obstruction of the primary 
flow by the secondary one. Lowering the injection temperature nearly broke even the amount 
of burning occurring. Overall, it appears that the slot injection of hydrogen at large angles 
to the primary flow is inadvisable relative to tangential or even low-angle injection. 

The performance of helium slot injection for tangential, 10°, and 45° injection angles is 
shown in Figure 103. For the cases tested, the use of helium as a coolant was always able 
to provide mitigation of ignition relative to the baseline. For both tangential and angled 
injection, helium was able to meet or outperform hydrogen when compared on the bases of 
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blowing rate and injectant temperature. This can be attributed to helium being an inert 
species. 

At high tangential blowing rates of helium, Figure 103(a) shows that helium was able 
to achieve similar performance to hydrogen at both coolant temperatures. Unlike hydro- 
gen, helium was not subject to the coolant burning observed in the hydrogen case. This 
allowed the helium coolant to experience a more gradual dropoff in performance compared 
to hydrogen as blowing was decreased. As the blowing parameter was decreased below 15%, 
low temperature helium experienced a slightly greater decrease in performance compared 
to high temperature helium. This was because the low temperature helium was traveling 
at a lower velocity versus the high temperature case. Therefore, more oxygen burning was 
occurring in the performance sampling region due to the shortened induction distances. 
It should also be noted that at low, 5% blowing levels, helium injection was still able to 
mitigate a small amount of autoignition, whereas hydrogen had excess burning at the same 
level of blowing. 

Over the range of blowing parameter values tested, 10° injection of helium shown in 
Figure 103(b) had worse performance relative to tangential helium injection. Compared to 
hydrogen, however, performance was better overall as there was no excess burning relative 
to the baseline case. Performance degradation as blowing was reduced was more gradual for 
low temperature helium versus the high temperature cases, which experienced a sharp drop 
in performance below B = 15%, with the general trend resembling a step function over the 
blowing parameter range tested. This was because at B = 15%, high temperature helium 
had a higher pressure of 12,700 Pa and achieved the “cushioning” effect initially desired for 
angled film cooling and mitigates ignition at the wall, shown in Figure 104(b). Lowering the 
blowing parameter below 15% eliminated this cushioning effect. Low temperature helium 
at B = 15%, with a lower pressure of 4600 Pa, did not have the same cushioning effect as 
shown in Figure 104(a). 

Injecting helium at a 45° angle to the main flow produced overall lowered performance 
relative to the 10° case, although there are some interesting features that can be observed in 
Figure 103(c). First, for low temperature helium, the performance at high levels of blowing 
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Tangential Helium Injection 


Te = 100 K 
Te = 740 K 



10 % = 20 % 15 % 10 % 5 % 

(a) 


10 deg Helium Injection 



45 deg Helium Injection 



(c) 


Figure 103: Changes in mass flow of oxygen burned for secondary helium injection is shown. 
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(b) 

Figure 104: Numerical shadowgraphs overlaid with water contours Yjj 2 o = 0.01, 0.02, ..., 
0.22 for 10 degree him cooling with helium. B = 15% and ambient coolant temperatures of 
(a) 100 K and (b) 740 K are shown. 

was similar to the same case for low angle injection. For high temperature helium, a well- 
defined local minimum in burning reduction was observed at lower blowing rates. This was 
due to the tradeoff between high amounts of shock-induced burning at higher blowing levels 
with lower levels of heat sink capability at low levels of blowing. Therefore, the decrease 
in oxygen burned between 20% and 10% blowing can be attributed to the decrease in how 
obstruction (and subsequent high temperature, shocked regions) and the increase in oxygen 
burned between 10% and 5% can be attributed to decreased cooling capacity. For the low 
temperature case, monotonically decreasing cooling capacity was witnessed over the entire 
range of blowing parameters because the coolant pressure at B = 20%, low temperature 
blowing (6620 Pa) was less than that for B = 10% high temperature blowing (8000 Pa). 
Thus, no local minimum was observed in the low temperature case. 

The performance of nitrogen is shown in Figure 105. Overall it was evident that nitrogen 
was a relatively poor performer in mitigating autoignition when compared to hydrogen and 
helium. Two drawbacks nitrogen had over hydrogen and helium at the same temperatures 
and blowing levels were its relatively low specihc heat with respect to mass and the low 
velocities with which it issued from the sonic slot. The first detriment impacted the ability 
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of the secondary fluid to act as a heat sink to the shear layer near the wall. The second 
detriment impacted the ability of the coolant to convect the hot flow downstream and delay 
the occurrence of ignition once the cooling ability of the injectant was saturated. 

Figure 105(a) summarizes the performance of tangential nitrogen injection. Overall, 
only modest mitigation of wall autoignition is achieved by nitrogen injection. Nitrogen does 
provide some small amount of heat sink capability, which is reflected in less oxygen burning 
with respect to the baseline. These gains in performance were roughly constant over the 
blowing parameter range with only marginal difference between low- and high-temperature 
nitrogen temperature. 

Injecting nitrogen at a 10° angle to the main flow resulted in worse performance relative 
to the tangential injection case as shown in Figure 105(b). None of the cases tested for this 
class of secondary injection were able to achieve over 10% reduction in the amount of oxygen 
burned relative to the baseline. Because of this, low-angled injection of nitrogen, while still 
providing marginal decreases in burned oxygen, was not a useful method of autoignition 
mitigation. 

Nitrogen injected at a 45° angle to the main flow resulted in a net increase in oxygen 
burning, although like for the previous cases, the increases in burning are marginal as shown 
in Figure 105(c). Therefore, it was evident that high-angle injection of nitrogen was also 
not a feasible method for autoignition mitigation, and furthermore increases the oxygen 
burning that was already occurring relative to the baseline case. 

3.5.3 Results of Wall Transpiration 

For wall transpiration, the cases indicated by checkmarks in Table 19 were simulated to 
convergence with the amount of oxygen burned analyzed using the same methodology as 
for the film cooling case. In addition, it was decided that any cases susceptible to the 
bleeding of main flow into the simulated plenum, a disastrous condition that would subject 
the internal transpiration system to hot gases, were failure cases. In addition, the plenum 
pressure was tracked by the simulation in order to arrive at an idea of the system demands 
for transpiration cooling. 
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10 deg Nitrogen Injection 
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45 deg Nitrogen Injection 



Figure 105: Changes in mass flow of oxygen burned for secondary nitrogen injection is 
shown. 
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Figure 106: Numerical shadowgraphs overlaid with water contours Yh 2 o = 0.01, 0.02, ..., 
0.22 for transpiration cooling with helium. Parameters B = 20% and W = 25% are shown. 

Of the cases shown in Table 19, only a subset of the helium cases were deemed to be 
feasible. Cases with hydrogen injection had the issue of causing autoignition near the leading 
edge of the transpiration patch behind the leading edge shock wave, shown in Figure 106. 
Injection of nitrogen as a coolant was demonstrated to be an ineffective cooling strategy 
in general for the same reasons of heat capacity as before, with the additional issue that 
the low plenum pressures associated with nitrogen tended to allow flow to bleed into the 
plenum. 

Of the helium cases tested, only the 25% width and 50% length cases were shown to be 
able to mitigate the presence of autoignition by any amount. However, only high values of 
blowing were shown to mitigate the presence of autoignition, and even then the gains were 
only minimal. Figure 108 shows that the burning characteristics of the best performing 
transpiration cases were similar to the baseline case without any cooling. This is because 
transpiration cooling injected fluid in the main stream at much smaller velocities relative to 
the main flow, severely reducing the cooling capacity of the fluid even for high specific heat 
species such as helium. Figure 107 shows the performance of helium injection for these cases. 
It was observed that a simple trend of increasing mitigation performance with increasing 
blowing parameter was achieved. At the highest blowing values, helium transpiration had 
performance only approaching that of tangential nitrogen film cooling. Because such poor 
performance was gained for the same amounts of mass flow as helium film cooling, it was 
concluded that transpiration cooling is not advantageous for autoignition mitigation for the 
conditions studied here. 
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Figure 107: Changes in mass flow of oxygen burned for secondary helium transpiration are 
shown. 
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Figure 108: Numerical shadowgraphs overlaid with water contours Yh 2 o = 0.01, 0.02, ..., 
0.22 for transpiration cooling with helium. B = 20% and patch sizings of (a) W = 25% and 
(b) L = 50% are shown. 
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3.5.4 Effect of Turbulent Schmidt Number on Autoignition 


It was observed in Chapter 2 that the turbulent Schmidt number affected the diffusivity 
of the strut injected fuel plume into the surrounding air. While this has an effect on the 
mixing efficiency achieved by the strut, combustion may also be affected due to the altered 
region of combustible fuel. In addition, the fuel may begin to enter the boundary layer 
sooner and, therefore, begin to burn at locations upstream of the baseline Set = 1.0 case. 
With these concerns in mind, it was necessary to evaluate the effect of the selection of Set 
on the predicted locations of autoignition and the efficacy of the implemented mitigation 
strategies. 

A subset of the cooling cases studied in the previous sections were selected in order 
to evaluate the effect of turbulent Schmidt number. The selected cases were taken from 
the tangential film cooling geometry and consisted of hydrogen and helium at the low 
(T e = 100 K) static temperature condition as these cases proved to be highest-performing 
of the Set = 1.0 cases investigated in the bulk of the study. All conditions and solution 
methodologies were the same with the only change being the change in Set to 0.25. 

Before evaluating the effect of turbulent Schmidt number on film cooling results, a 
reacting baseline case at the lower turbulent Schmidt number was computed without any 
film cooling. This is because the results here are presented in terms of a relative change 
in oxygen burned. For a fair comparison, the change in oxygen mass burned should be 
relative to a baseline case run with the same turbulent Schmidt number. After re-running 
the reacting strut case (without cooling) with Set = 0.25, the oxygen mass burned changed 
from 5.679e-3 kg/s to 8.728e-3 kg/s, a 54% increase relative to the Set = 1.0 case. 

The changes in mass of oxygen burned for each Schmidt number for the selected hydrogen 
and helium cases are shown in Figure 110. The first comparison case was low temperature 
hydrogen injection. It is observed that lowering the turbulent Schmidt number overall 
increased the percent change in oxygen burned, but this figure is relative to the respective 
baselines. Figure 109 shows the actual burned mass flow rate of oxygen for each of the 
tested cases. It is observed that at high secondary hydrogen mass flow rates, lowering the 
turbulent Schmidt number slightly decreased the amount of oxygen burned while at low 
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mass flow rates the amount of oxygen burned increased. For secondary helium injection, 
lowering the turbulent Schmidt number increased the oxygen burned across all blowing 
rates, although at 20% the values are very close. Overall, it is concluded that changing the 
turbulent Schmidt number did have an effect on the burning characteristics observed, but 
did not significantly impact the high performance of tangential film cooling of hydrogen and 
helium. 

3.5.5 Effect of Chemical Model on Autoignition 

It is useful to understand the sensitivity of the numerical simulation to the selection of 
chemical model in order to increase confidence in the results obtained in this study. The 
strut and stinger baseline autoignition cases were rerun with the O Conaire [92] chemistry 
model with and without reactions involving H 2 O 2 . This gives knowledge not only to the 
sensitivity of the problem to a different, modern mechanism, but also to the importance of 
hydrogen peroxide to the reacting system. 

Contours of water production and mass fraction downstream of the strut and stinger 
injectors are shown in Figures 111 and 112, respectively. The overall morphology of the 
contours agree with those from the baseline simulations with Jachimowski chemistry. The 
oxygen consumption for the strut case decreases from 0.005728 kg/s with the Jachimowski 
baseline to 0.005329 kg/s, a 7.2% decrease drop relative to the mean. Adding hydrogen per- 
oxide chemistry from the O Conaire model affects a 6.7% decrease relative to the baseline. 
For the stinger case, the oxygen consumption decreases from 0.0086 kg/s with the Jachi- 
mowski baseline to 0.008415 kg/s, a 2.17% decrease drop relative to the mean. Adding 
hydrogen peroxide chemistry from the O Conaire model affects a 2.12% decrease relative to 
the baseline. These results indicate that while the solution is sensitive to a small degree to 
the selection of chemical model, it is not very sensitive to the presence of hydrogen peroxide 
chemistry at the simulated conditions. Therefore, it may not be worth solving for hydro- 
gen peroxide chemistry depending on the additional number of equations (and computation 
time) that will be required. 
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Tangential Helium Injection 
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Figure 109: Changes in mass flow of oxygen burned for tangential, secondary hydrogen 
and helium injection is shown. 
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Tangential Hydrogen Injection 




(b) 

Figure 110: Percent changes in mass flow of oxygen burned for tangential, secondary 
hydrogen and helium injection is shown. 
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(b) 


Figure 111: Stmt injector contours of (a) water production and (b) water mass fraction 
in the flow domain of interest using the O Conaire mechanism without hydrogen peroxide. 
Flow is from the upper left to the lower right. 
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(b) 


Figure 112: Stinger injector contours of (a) water production and (b) water mass fraction 
in the flow domain of interest using the O Conaire mechanism without hydrogen peroxide. 
Flow is from the upper left to the lower right. 
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3.6 Chapter Summary 


The current study of the mitigation of autoignition on the forebody of a hypervelocity vehicle 
focused on the use of active wall cooling. Two methods of wall cooling were considered: slot 
injection (both parallel and angled) and wall transpiration. Both strategies were applied 
as boundary conditions, although the transpiration modeling used in the study allowed 
for physical reproductions of the presence of suction and blowing on the same region of 
transpiration. For each method of cooling, three species were tested: hydrogen, helium, and 
nitrogen. Hydrogen was chosen because of its systems benefits (by using a single fueling 
system) and excellent heat capacity at the cost of being a combustible coolant. Helium was 
chosen as an alternative because of its marginally lower heat capacity versus hydrogen with 
the additional benefit of being an inert species. Nitrogen was chosen as a third alternative 
because it is also inert and has improved storage density over helium, although its heat 
capacity is very low compared to the other two alternatives. Both film and transpiration 
cooling strategies were evaluated at various levels of blowing. 

Baseline reacting cases of the stinger, unswept ramp, and strut injectors from Chapter 
2 were simulated without any wall cooling in order to provide a baseline case for compari- 
son. The study of autoignition mitigation focused on the strut injector because transverse 
injection through a stinger proved itself to be susceptible to extensive amounts of ignition 
originating away from the wall while encapsulating the fuel plume. The unswept ramp 
injector also exhibited redundant ignition characteristics at the wall compared to the strut, 
just further down stream. Therefore, it was found that the strut provided an excellent 
geometry for study, with any solution able to potentially be applied to ramp injection. 

Film cooling was tested first for hydrogen, helium, and nitrogen over blowing parameters 
ranging from 5% to 20%. Two exit temperatures of 100 and 740 K were tested to determine 
the effect of fluid total temperature on the cooling characteristics. It was found that at 
high blowing rates and low temperatures that hydrogen was able to eliminate all of the 
wall autoignition occurring for the strut case. This effect decreased in performance both as 
coolant temperature increased and/or blowing rate was lowered. Performance also degraded 
as angled injection was introduced as the coolant itself became susceptible to burning. In 
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general, helium was able to achieve similar levels of performance as hydrogen and sustained 
better performance at lower levels of blowing due to its inert character. Overall, nitrogen was 
a poor choice in mitigating autoignition due to its diminished heat capacity with respect 
to mass. Therefore, it was concluded that tangential hydrogen and helium cooling are 
effective methods for autoignition mitigation. If high cooling mass flow rates of hydrogen 
are required to mitigate autoignition, fueling splits between the primary and film cooling 
injectors may be considered to target the optimum system global equivalence ratio. 
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CHAPTER IV 


TWO-DIMENSIONAL COMBUSTOR CALCULATIONS 


The intentional shocking of the premixed fuel-air mixture from the forebody flow is the 
primary mechanism by which ignition occurs in a shcramjet concept. While wall burning 
can contribute to the consumption of fuel in the combustor, it is the shockwave that traverses 
the entire flow that produces the bulk rise in temperature ideally resulting in nearly uniform 
combustion. Typically the shock wave is generated by the cowl lip due to the flow turning 
at that location. In some concepts, however, the shock wave that raises the temperature of 
the flow may be generated by a wedge placed at the wall of the internal flowpath. 

Depending on the point on the trajectory the vehicle is located, the separation distance 
between the shock and reaction fronts will vary to a large degree. Higher Mach numbers 
will yield correspondingly higher post-shock temperatures, thereby decreasing the ignition 
delay time experienced by the shocked mixture. Therefore, as the Mach number of the 
trajectory is increased, the shock and reaction fronts will become more closely spaced, 
eventually merging. As the two fronts approach each other, there is the possibility of 
unstable interaction with each other as the normal velocity would be subsonic, but this has 
not yet been observed in simulation. At low enough Mach numbers, there will either be 
little or no ignition or the ignition delay time will be long enough to neutralize the benefits 
of the shcramjet due to the requirement of a long combustor. These Mach numbers will 
functionally define the lower limit of feasibility for a shcramjet. 

In this chapter, attention will be given to the simulation of premixed, shock-induced 
combustion in a notional combustor of arbitrary length. For this study, it was assumed that 
the cowl wall provided the necessary turning of the flow in order to produce a shock wave 
that raised the temperature of the gas above its autoignition limit. Here the angle of turning 
was limited to roughly 11 degrees, a typical used total turning for such flows. This study is 
important in establishing quantitative performance metrics and qualitative characteristics 
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of the combusting flow. Information gained from this study included the determination of 
optimum reaction length as a function of freestream conditions on a notional trajectory as 
well as stability characteristics along the trajectory in the vicinity of the cowl leading edge. 

4-1 Modeling and Simulation 

4.1.1 Inflow Profiles 

Evaluation of burning characteristics in the notional combustor used two types of profiles 
as the inflow: uniformly premixed fuel-air and profiles derived from the mixing simulations 
of Chapter 2. For the cases derived from the mixing simulations, the outflow plane of 
the simulation was saved and used as the inflow for the notional combustor simulations. 
This allowed the combustion characteristics of selected cases to be evaluated relatively 
inexpensively without needing to exhaustively simulate a notional combustor for each mixing 
case. The fuel-air profile of the strut from the mixing study is shown in Figure 65 of Chapter 
2 for x = 0.713 m. Note that nonreacting profiles were selected, so the implicit assumption 
is that no reaction had occurred until cowl closure, which was shown to be not true in 
Chapter 3. This assumption was made as a simplifying measure. 

For studying uniformly premixed combustion, integer freestream Mach numbers from 
a Mach 10 to Mach 20, constant 71,820 Pa (1500 psf) dynamic pressure trajectory were 
simulated. Inflow conditions were taken from a custom-made trajectory code that takes the 
oncoming flow conditions at a particular trajectory point, shocks it through two turning 
angles according to the notional forebody used, and then ideally mixes the airflow with the 
fuel while conserving mass, momentum, and energy. More details about the methodology 
are described in Appendix C. At each trajectory point the fuel flow supplied to the forebody 
air was varied in order to provide the required equivalence ratio at that point. 

4.1.2 Computational Domain and Methodology 

The topology used in this study consisted of a short runway for grid transition from the 
outflow plane of the mixing studies to one that was resolved on both walls of the internal 
flow path. The lower wall of the topology consisted of an isothermal wall through the 
entire length. The upper boundary had an extrapolation boundary condition until the 
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Figure 113: The grid used in the combustor study is shown with every fourth grid point 
plotted normal to the walls. 

point of cowl closure, where an isothermal wall was defined and continued to the end of the 
computational domain. The sides of the computational domain are defined as symmetry 
boundary conditions. The cowl wall was angled at 10.8° to the oncoming flow in order 
to match the total turning experienced by the notional forebody. Figure 113 shows the 
topology used in this study. The number of grid points in the topology totaled 1.46 million 
with 1.41 million points located behind the cowl closure location. 

Past the shoulder of the topology was a constant area duct that served as a notional 
combustor. In an actual flight vehicle the combustor walls would gradually diverge in 
order to combat thermal choking, but consideration is not given to that in this study. 
The length of the combustor was arbitrarily defined to be one meter to allow combustion 
to substantially complete at lower Mach numbers. The total length was chosen arbitrary 
because the purpose of this study was not to evaluate the performance of an actual design 
combustor, but rather to comparatively evaluate the performance characteristics of different 
inflow mixing profiles. 

For simulations of premixed flow stability in the locality of the cowl leading edge, a 
refined grid was constructed using the same geometry as the full simulation, but cut off a 
short distance past the shoulder of the geometry. The grid had dimensions 256 x 587 points 
axially and transversely and is simulated in two dimensions. The topology is shown in 
Figure 114. The walls, as with the full simulation, were set to be at a constant temperature 
of 500 K. 

The flow was solved using the computational methodology outlined in Section 2.1, with 
the chemistry using Jachimowski’s reduced model described in Section 3.1.1. 
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Figure 114: The grid used in the stability study is shown. 

4.1.3 Performance Metrics 

Combustion efficiency and thrust potential were the quantitative metrics used to track the 
performance characteristics of different inflow conditions. Combustion efficiency was used 
as a way to track the progress of the combustion within the notional combustion chamber 
by tracking the consumption of the species not in excess. As discussed in Section 3.2.2, 
the definition of combustion efficiency varies. Here it was defined as the consumption of 
oxygen mass relative to the inflow oxygen mass. The stream thrust potential, used as a 
loss tracking metric in the previous study, was used here to track both the thrust benefits 
gained by combustion as well as the viscous and shock losses throughout the flowpath. This 
use case further supports the use of thrust potential throughout the flowpath as it gives 
an honest assessment of the balance between the beneficial and detrimental mechanisms 
occurring within a flow. 

Using thrust potential allows the determination of optimum flowpath lengths for com- 
ponents such as engine modules [69] because the general thrust potential trend peaks at the 
point where gains due to combustion begin to lose to viscous losses near the walls of the 
flowpath. Figure 115 shows a notional plot of thrust potential of gas burning in a combustor 
as a function of distance along the flowpath. Because stream thrust potential is a measure 
of the stream thrust obtainable by isentropically expanding a particular cross-section to 
a reference nozzle exit, the optimum point on the plot represents the ideal point for geo- 
metric expansion to occur. Therefore, the distance between the beginning of the notional 
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Figure 115: A notional trend of thrust potential versus distance in a combustor. 


combustor and the local maximum is also referred to as the optimum combustor length. 
This length is related to the ignition delay time of the fuel-air mixture. Therefore, higher 
Mach numbers will demand a shorter combustor length because burning occurs relatively 
rapidly due to the shortened ignition delay times. 

4- 2 Results and Discussion 

Integer points between Mach 10 and Mach 20 on a constant 71,820 Pa (1500 psf) trajectory 
were simulated over a one meter length constant-area duct after being shocked via turning 
by a notional cowl. The premixed conditions were obtained by shocking the freestream 
through two turns and then mixing the air with fuel at stoichiometric proportions. Figure 
116 shows contours of water mass fraction overlaid with lines of constant pressure for Mach 
numbers 11, 15, and 20. It was observed that as freestream velocity increased, the separation 
distance between the shock and reaction fronts decreased accordingly. This was due to the 
fact that larger upstream Mach numbers induce a larger temperature rise which had the 
effect of decreasing the induction distance through decreased ignition delay times (even 
though velocity increases, which would have the effect of delaying ignition for different 
cases having the same temperature). 

The Mach 11 solution did not exhibit significant shock-induced combustion, but rather 
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(c) 

Figure 116: Contours of water mass fraction overlaid with lines of constant pressure for 
Mach numbers (a) 11, (b) 15, and (c) 20. 
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the high temperatures of the boundary layers caused the onset of combustion. This is 
because the temperature immediately downstream of the oblique shock wave was roughly 
820 K, too low for autoignition. The combustion zones at the wall are increased in size due 
to the presence of expansion waves in the flow reflecting off of the walls of the flowpath, 
with both combustion zones eventually merging downstream. At higher Mach numbers, 
for example Mach 15 in Figure 116(b), a shock-induced reaction front formed due a higher 
temperature behind the shock wave, approximately 1180 K in this case. As freestream Mach 
number was increased, the shock and reaction fronts eventually began to merge. Note also 
that as freestream Mach number was increased, the mass fraction of water in the burned 
flowfield is lessened mainly due to the very high temperatures in the flow. For example, the 
post shock temperature of the Mach 20 case is 2200 K, causing the dissociation of water. 

Decreasing the distance between shock and reaction fronts generally had a beneficial 
effect on the overall system because it allowed for shorter optimum combustion lengths. 
This is because the fuel-air mixture ignited and proceeded toward completion in shorter 
distances. However, combustion efficiency is not an ideal metric to use for designating the 
optimum combustion length; stream thrust potential is a better metric for such purposes. 
Figure 117 shows plots of thrust potential and combustion efficiency for the Mach 12 case. 
The combustion efficiency curve is characterized by a shallow initial rise followed by a rapid 
rise as oxygen is rapidly consumed. Eventually the rise in combustion efficiency tapers off 
past the knee in the curve. Now examine the stream thrust potential curve, which increases 
in value due to combustion and reaches a maximum at the critical point where viscous losses 
begin to overcome the gains from further combustion. 

When given only combustion efficiency, an analyst could pursue a few options for defining 
an “optimum” combustion length. First, he may define a “threshold” combustion efficiency, 
say 80%. This may not be a good idea for multiple reasons. The first is that it is unknown 
if this value is attainable for a given condition due to the high temperatures experienced in 
hypervelocity flight. Second, if this value is attainable, it may demand unreasonably long 
combustor lengths. This can be observed in Figure 117, which shows the combustion curve’s 
knee at a combustion efficiency of 0.70, past which point combustion efficiency increases 
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X [m] 

Figure 117: Thrust potential and combustion efficiency for the Mach 12 case. The shoulder 
of the geometry is where x = 0. 
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slowly. Either way, the optimum point in the thrust potential is observed to coincide with 
a combustion efficiency value of 0.76, so designing the combustor past this point produces 
increasing losses anyways. Another option for defining combustor length by using the point 
of the knee in the curve is also not ideal because it is observed from the thrust potential 
curve that greater gains can be earned by a greater length. 

Figure 118 plots the optimum combustion length for different points in the notional 
71,820 Pa (1500 psf) trajectory. The curve of optimum combustion length can be split 
into three regions. The first, below Mach 13, is dominated by boundary layer burning and 
does not exhibit a rapid drop in optimum combustor length as a function of freestream 
Mach number. Between Mach 13 and Mach 17, a rapid drop in optimum combustion length 
occurs as a shock-induced combustion front travels upstream toward the cowl shock. Above 
Mach 17, the trend stops decreasing as the combustion front dwells in the vicinity of the 
initial shock wave. The combustion efficiency at the location of maximum stream thrust 
potential is also plotted and shows that in general as freestream Mach number is increased, 
the combustion efficiency at the optimum point decreases. This confirms the conclusion 
that combustion efficiency is not an ideal metric to use for finding the optimum combustion 
length as the “optimum” combustion efficiency changes as a function of freestream Mach 
number. 

The trend in optimum combustor length confirms that as freestream Mach number 
increases, the distance required in order to extract maximum stream thrust potential de- 
creases. For the particular notional trajectory considered, there were a wide range of opti- 
mum combustion lengths. At Mach 11, which also corresponds to the low end of feasibility 
for the conditions tested, the optimum combustion length is roughly 0.8 m, while at Mach 
20 the length shortens significantly to 0.1 m. During the design of an actual vehicle flow 
path with fixed geometry, a plot like this may be used to determine the minimum required 
combustion length by identifying the point on the curve where shcramjet concept is to begin 
operation. On the other hand, if a variable internal geometry is used, a plot such as this can 
provide a schedule for when to begin the expansion of the flow to the nozzle exit in order 
to extract the largest possible performance from the flow. The large variation in required 
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Figure 118: Optimum combustor length and combustion efficiency at that point at different 
trajectory points. 

combustor length suggests that variable internal geometry will be an enabling technology 
for actual flight shcramjets. 

The flowfield input with a nonreacting strut injector cowl closure profile is shown in 
Figure 119. Three longitudinal slices are shown at the centerline plane, the symmetry plane 
between injectors, and the plane positioned between the two symmetry planes. Here it was 
observed that the shock-induced combustion occurs relatively rapidly compared to the Mach 
12 premixed case. This was primarily due to the elevated temperatures in the strut flowfield 
facilitating shorter ignition delay times. It was also observed that the reaction front forms 
at the location of the shock reflecting away from the shoulder region of the geometry. Upon 
examination of the stream thrust potential and combustion efficiency trends, the stream 
thrust potential was observed to reach its peak at 0.681 nr, a bit upstream of the premixed 
case, with a combustion efficiency of 0.78. 

In order to determine the presence of flow instabilities, integer freestream values on the 
Mach 15 to 20 trajectory were simulated in the vicinity of the cowl closure location with a 
refined grid. This was done because combustion instabilities occurring in a flow are sensitive 
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(c) 

Figure 119: Filled contours of water mass fraction superimposed with lines of constant 
pressure for the Mach 12 strut flowfield at y = (a) 0, (b) 0.0172, and (c) 0.0345 m 
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to grid resolution and simply simulating the entire notional combustor flowpath at a high 
grid density would not be practical given where the shock-reaction front instabilities would 
occur (i.e. near the shock wave). 

Figure 120 shows the shock-reaction front system for Mach numbers 15 through 20 for 
the more refined simulation in the vicinity of the cowl leading edge. For the cases shown, as 
the reaction front approaches the shock front no instabilities are visible in the shock layer 
or in the reaction zone. Some shock instabilities have been reported in the literature at 
higher wedge angles, but at the low turning angles characteristic of airbreathing propulsion 
there may be no instabilities for engineering purposes. 

4-3 Chapter Summary 

The purpose of this chapter was to examine the properties of the reacting flow system in an 
internal flowpath representing a notional combustor. Both ideally premixed fuel-air as well 
as the strut fuel-air profile as extracted from the outflow plane of the mixing study were 
evaluated. Two grids were used for the study, one that simulated the cowl closure shockwave 
as well as a constant-area duct with a length set to 1 m as well as a more refined grid in 
the vicinity of the cowl closure region in order to ascertain the presence of instabilities in 
the shock-reaction front region. 

In simulating a constant area duct with reactions, the point of maximum stream thrust 
potential was sought. This maximum occurs at the location where gains due to the burning 
of the fuel and air were balanced by the losses due to viscous losses at the wall. Such a 
point represents the location where the flow would ideally begin to expand to the nozzle. It 
was noted, however, that this optimum location varied significantly depending on the point 
on the trajectory the vehicle is currently located. This underscores the importance of the 
development of variable internal geometry in order to accommodate the maximum possible 
performance from the engine. 

In order to evaluate the presence of instability for cases where the reaction front and 
shock were closely spaced, Mach numbers between 15 and 20 were simulated in the vicinity 
of the notional cowl leading edge with a fine grid. While the contour of the formation 
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Figure 120: Filled contours of water mass fraction superimposed with lines of constant 
pressure for Mach numbers (a) 15, (a) 17, and (c) 20. 
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of H2O closely followed the cowl shock wave at high Mach numbers, no instabilities were 
observed in the flow. This may be due to the low turning angle of the flow relative to other 
cases in the literature that show instabilities at higher turning angles. 
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CHAPTER V 


QUASI-ID SHOCK & REACTION FRONT STABILITY 


Instability in shock-induced combustion has been a subject of experimental attention since 
the early 1960’s. Early classic experiments were conducted by Ruegg and Dorsey [96], 
Hartmuth F. Lehr [10], McVey and Toong [11], and Alpert and Toong [12] of PMSIC ahead 
of spherically blunted projectiles. These experiments identified three important stability 
regimes for such flows: stable; regular, unsteady; and large-disturbance, unsteady. 

Stable Regime The stable regime is characterized by a reaction front and shock front 
with standoff distances that are invariant with time. Figure 121 shows a notional illustration 
of the reaction and shock front locations on the stagnation streamline as a function of time. 
Because combustion in this flow regime is not time varying, it is the most desirable regime 
for shcramjet combustion. 

Regular, Unsteady Regime The regular, unsteady regime is characterized by peri- 
odic, coupled oscillations in both the reaction and shock fronts. Figure 122 shows a notional 
illustration of the reaction and shock front locations on the stagnation streamline as a func- 
tion of time. The currently-accepted mechanism for this behavior is as follows. First, a 
upstream-traveling compression wave interacts with the bow shock by pushing the shock 
forward and reflecting a contact discontinuity downstream. As the induction time on the 
upstream side of the contact discontinuity is lower due to higher temperatures compared to 
the downstream side, a new reaction zone forms upstream of the old reaction front. This 
causes compression waves to be generated and travel downstream toward the body and 
upstream toward the shock, eventually generating another contact discontinuity. Once the 
two reaction zones merge, an expansion wave travels upstream and relaxes the shock wave 
back downstream. This cycle repeats itself indefinitely and at regular intervals. 
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Figure 121: Reaction and shock front locations on the stagnation streamline as a function 
of time for the stable regime. 
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Figure 122: Reaction and shock front locations on the stagnation streamline as a function 
of time for the regular, unsteady regime. 
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Figure 123: Reaction and shock front locations on the stagnation streamline as a function 
of time for the large-disturbance, unsteady regime. 


Large-Disturbance, Unsteady Regime The large-disturbance, unsteady regime is 
characterized by large-amplitude disturbances in the reaction zone that cause it to merge 
with the shock front periodically. This regime has been studied in the past by Matsuo, et 
al. [105]. Figure 123 shows a notional illustration of the reaction and shock front locations 
on the stagnation streamline as a function of time. In this mechanism, a deflagration-to- 
detonation transition (DDT) [106, 107] occurs at the point where the downstream-traveling 
contact discontinuity interacts with the reaction front. At this point, a retonation wave 
travels downstream toward the body while a detonation wave travels upstream and merges 
with the shock wave. Eventually the reaction and shock fronts separate and the process 
repeats itself. This behavior is very large amplitude and in some cases has been observed 
to be non-periodic. Due to the violently unstable nature of this flow regime, it is highly 
undesirable for a shcramjet flowpath. 

Past numerical studies of PMSIC ahead of a blunted projectile [108, 105, 109, 110, 
111, 112] have been carried out on two dimensional, axisymmetric grids, shown in Figure 
124. Because the PMSIC phenomenon requires a reacting CFD simulation with typically 
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Figure 124: Schematic of the current state of the art for solving flow ahead of a blunted 
projectile. 

high grid density, these simulations can take a very large amount of computational work 
to arrive at a converged solution. Because of the time-varying nature of such flows, a 
solution is considered to be converged when the initial, transient behaviors of the flow 
have passed. Due to the long clock time required for convergence and the time required 
at the converged condition in order to get good data for characterizing the flow, it can be 
difficult to understand the behavioral trends of this flow as a response to adjustments in 
the freestream condition (e.g. equivalence ratio, velocity, density, etc.). 

An advancement to the state of the art is proposed by using a quasi- ID formulation 
for PMSIC ahead of a blunted projectile. Previous experimental and numerical work has 
observed that the unsteady features of such flows are symmetric [96], suggesting that the 
driving mechanism is emanating around the stagnation streamline. If the flow solution is 
carried out only along the stagnation streamline, the computational costs associated with 
the multi-dimensional simulations can be greatly reduced with a minimal price paid in terms 
of accuracy. A schematic of such a solution domain is shown in Figure 125. 

In the quasi- ID formulation, a local similarity postulation is applied to the Navier- 
Stokes equations. The results is a system of equations that are a function in space only 
as distance from the surface of the projectile. The term quasi-ID refers to the fact that 
the simplified governing equations account for the downstream influence by including the 
gradient in the tangential velocity. With the quasi- ID simulation framework in place, it is 
possible to rapidly determine the stability characteristics of the flow in response to variation 
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Figure 125: Schematic of a new way for solving flow ahead of a blunted projectile. 


in inflow equivalence ratio, velocity, density, fuel type, and vitiation. Variation of inflow 
equivalence ratio as a function of time may also be simulated in order to model the effect 
of vortices. 

The quasi- ID formulation has high value for fundamental PMSIC research in quickly 
identifying the specific mechanisms of different regimes of stability with greater detail than 
what may be practical with multidimensional simulations. New trends can also be quickly 
identified due to the ability to run many simulations in a short span of time. Large design 
spaces can be populated with defined regions where the different stability mechanisms exist. 
Beyond fundamental research, such a capability can be useful for practical engine design. 
For example, corridors (as a function of velocity and density) can be identified where blunt 
body features in the flow will not generate the instabilities seen in previous experimental 
and numerical work. 

This purpose of this chapter will be to summarize the development, validation, and 
application of the quasi-ID form of the Navier Stokes equations in order to enable the 
rapid analysis of shock-induced combustion due to premixed fuel-air in blunt-body flow. 
After deriving the reduced equations, a description will follow of the physical modeling that 
supports the solution of the equations in a custom-built analysis tool. The tool was then 
validated against theory and experiment. Validating the reduced framework was a primary 
deliverable of the study. The framework was also be applied to conduct parameter studies 
previously unfeasible due to the computational demands of axisymmetric solutions. 
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Figure 126: The flow schematic and coordinate system for the quasi-ID equations are 
shown. 


5.1 Derivation of the Quasi- ID Form of Equations 

5.1.1 Governing Equations 


The axisymmetric Navier-Stokes equations and species conservation equations in the body- 
oriented coordinate system shown in Figure 126 can be written as 


<9U <9M <9N ^ 

aT + ^7 + + Q 


( 12 ) 


The state vector (U), flux vectors (M,N), and the source term vector (Q) are expressed in 
non-dimensional form as 
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where A = 1 + nk and /3 = r + n cos a. In these equations, variables have been non- 
dimensionalized with respect to free stream velocity, temperature, density, and body radius. 

In order to arrive at the quasi- ID form of equations, a local similarity postulation was 
used in the stagnation region. Under this approach, flow variables were expanded in a series 
about the axis of symmetry (see Appendix A), resulting in the state equations 
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Note that the retention of the second order term in the pressure was warranted by the 
order of magnitude analysis. The expansions given in Equation 17 have had a small-angle 
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approximation applied due to the result that s = a for a non-dimensionalized sphere. The 
solution variables in the first truncation have a simplified notation where terms in the 
expansion equal to zero have been removed and the remaining solution variables and their 
derivatives have simply been given the names of the expanded state variable with subscripts. 
Here u± is actually the first derivative of velocity with respect to s and P 2 is the second 
derivative of pressure with respect to s. All other subscripted values in the first truncation 
are equal to their respective state variables. 

Retaining terms only to the first truncation shown in Equation 17 and substituting 
into Equations 13 - 16, the following set of equations is obtained for a sphere by equating 
the coefficients of like powers of a. Note that on the stagnation streamline a = 0. The 
governing equations are reduced to the form 
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where now the solution is a function only of time and distance normal to the surface. Note 
that the primes have been added to indicate the quasi- ID form of the vectors. The state 
vector has the form 
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The flux vector is split into contributing terms by 


N' = N- + N' v + N; 


(20) 


where each term represents the inviscid, viscous, and chemical contributions to the flux, 
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respectively. These terms are given by 
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An additional equation is required for the solution of p 2 , given by 
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which comes from the s - momentum equation. This equation is decoupled from the Navier- 
Stokes equations and can be solved separately in a non-lagging manner. 

The heat flux in the solution is composed of viscous and diffusive components, which 
contribute to the total heat flux according to 


Q2 — Q2,visc T Q2,diff (28) 

where the nondimensional viscous and diffusion heat fluxes have the form 

Hi dh\ 
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An additional step can be taken to convert the equations in the physical domain to the 
computational domain in order to streamline the handling of nonuniform grid spacing or 
grids that change in time. More details on this transformation can be found in Appendix 
B. 

5.1.2 Boundary and Initial Conditions 

The numerical solution of this problem requires the application of specific boundary con- 
ditions at the free stream wall as well as a sensible initial condition across the solution 
domain. The free stream boundary conditions are given by 



tti(oo) = pi(oo) = Ti(oo) = 1 
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with the boundary conditions for pressure and enthalpy given by the enthalpy and perfect 
gas equations. The wall conditions are set to be equal to those in the neighboring grid point 
with the exception of velocity. For the inviscid case v\ is set to zero at the wall and for the 
viscous case both u\ and v\ are set to zero. 

The initial condition across the solution domain should not contain any sharp gradients 
in order for the solution to be convergent. Therefore, the initial condition of each variable 
is assumed to be monotonically increasing or decreasing depending on the initial guess for 
the wall conditions. 

5. 2 Physical Modeling 

This subsection describes the thermochemical underpinnings of the tool developed to study 
the quasi-ID problem, the Quasi-ID solVer with Reactions (QUIVR). Although the solu- 
tion domain of the problem is very simple and described elegantly by a modified form of 
the Navier Stokes equations, computing the thermochemical properties of each cell requires 
the implementation of models describing reaction rates, equilibrium constants, viscosity, 
conductivity, specific heat, enthalpy, entropy, and Gibbs free energy. While it is possible to 
leverage a computational framework to compute these quantities from first principles (e.g. 
statistical mechanics), it is much more practical to implement models using polynomial 
curve fits. Furthermore, while open source packages (i.e. CEA and Cantera) may be called 
upon to compute these properties, enhanced portability and lower overhead is gained by 
programming these models specifically for the tool. 

5.2.1 Nonequilibrium Chemistry 

The most general description of the chemical evolution of the flow is that it is in chemical 
nonequilibrium. Time-accurate changes in chemical composition are addressed through 
the Navier-Stokes equation that conserves mass, with one equation for each species. The 
production rate for species i appearing in Equation 26 is given by the formula 

f)c • 
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here with the left hand side having the units kg/m 3 /s. The concentration time rate of 


change of species i is given by 
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where the right hand side is the series sum of the gain and loss terms of species i for reactions 
1, N. Casting the time rate of change of species concentration in terms of gain and loss 
terms is useful for the a-QSS solution procedure outlined in Section 5.3.3. The gain and 
loss terms have the form 


Gij = v'/jqj (33) 

Lij = v ijQj (34) 


with the rate-of-progress variable for reaction j having the form 


% = k f,i n - kb J n ^ ( 3s ) 

i i 

for reactions not containing a third body reaction and 

% = k fj ^ (n c 7 j - k fj n ^ ( 3e ) 

for reactions using a third body. Mechanism models using third body efficiencies typically 
list only species whose third body efficiencies are greater than unity for a particular reaction. 
All remaining unlisted species were assumed to have a third body efficiency of unity. Note 
that some mechanism models (e.g. Baulch, et al. [113]) use reactions that employ Troe 
parameters in order to introduce pressure dependency to certain reactions. This additional 
level of complexity was not accounted for in the chemistry implementation of the quasi- ID 
tool and was left for future work. Until then, chemistry models using Troe parameters are 
incompatible with the current implementation. 

The forward reaction rate coefficient for each reaction was assumed to be of Arrhenius form, 
given by 

kf = C,T°< exp (0) (37) 

where the reaction rate constants are generally listed in any given mechanism model along 
with the corresponding reaction. The backward rate coefficient is obtained simply by the 
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relation 


kb — 


ht 

K c 


(38) 


where the equilibrium constant in terms of concentration is calculated at a given tempera- 


ture from thermodynamic data and given analytically by the equation 


K c 



-An,' 

exp 



(39) 


where 


N N 


AG T - Y ~ Y 

(40) 

2=1 2=1 


N N 


A,, .i = Y , - Y "h 
2=1 2=1 

(41) 


The presence of the reference pressure, p°, in Equation 39 ensures that the units are con- 
sistent between the coefficient term and the exponential term. 


5.2.2 Viscosity and Thermal Conductivity 


The viscosity and thermal conductivity of the gas mixture was not computed by a using 
a polynomial form (although such a form does exist [114]). Rather, it has been proposed 
that the best compromise between transport model complexity and accuracy is given by the 
forthcoming equations [95] . The dynamic viscosity and thermal conductivity of a monatomic 
species in a pure gas are given by 

26.693 y/rriiT 
,U = afnt 2 ’ 2 )* 

15 R u / . Q, 

ki = m (43) 

4 mi 

The dynamic viscosity for a molecule is given by the same equation as for a monatom while 
the right hand side of the molecular thermal conductivity is given by 




m 





(44) 


where fi n t = 1.32, viscosity is in units of micro-Poise, and thermal conductivity is in units 
of cal/s/cm/K. The value of the Lennard-Jones (12-6) 1 collision diameter that appears in 


twelfth-power repulsion; sixth-power attraction 
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Equation 42 is species-dependent and is tabulated by Svehla [115]. The collision integral 
appearing in the same equation is underpinned by the bimolecular collision physics of the 
system and is given by 

q( 2 ’ 2 )* = ri( 1 ’ 1 )*A* (45) 

where 

= ao + ai T* + a 2 (T*) 2 + a 3 (T*) 3 + a 4 (T*) 4 (46) 

A* = ao + aiT* + a 2 (T*) 2 (47) 


The coefficients in in Equations 46 and 47 are chosen depending on the range of values in 
which the reduced temperature is located and are summarized in Table 20. The coefficients 
are also independent of the species being considered; however, the reduced temperature in 
these equations is species-dependant and is given by the relation 


rj~i^ 


e/k B 


(48) 


where e/k B is the species-dependant Lennard- Jones (12-6) energy of attraction divided by 
the Boltzmann constant. 

For a multi-species mixture, the mixture values of dynamic viscosity and thermal con- 
ductivity must be calculated. Mixture viscosity is given by Wilke [116] as 

N 


E Ah 
i= 1 1 + Xi ^3= 1 


Ah 

=1 ,j^i Xj^ij 


(49) 


where 


4 


l^ij — 


1+ (Mi) 2 ( Hh ' 4 
' ' Hj ) \m,i 


l | rrij 


(50) 


The mixture thermal conductivity is approximated by Coffee and Hirnerl [117] as 


km.ix — „ 


-r 


E*A + (Ef 


V 7=1 


(51) 
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5.2.3 Diffusion Velocity 


The calculation of the species diffusion velocity arises in the species conservation equation. 
The general method for determining diffusion velocities requires the solution of systems 
of equations; however, a simpler method is commonly used in practice. An approximate 
calculation based on Fick’s law of the first N — 1 species velocities is given by the equation 


Vi 



Vfi 


(52) 


while the N th species diffusion velocity — typically computed for the species in excess — is 
given from conservation form as [118] 


N - 1 


vn = — 


fN 




(53) 


%— 1 


which ensures that the species diffusion fluxes sum to zero. 

The effective mixture binary diffusion coefficient of species i is given by 


Dim. — 


1 -Xt 
xj 

Dij 


(54) 


which combines the species binary diffusion coefficients. The binary diffusion coefficients 
are defined by collision integrals that must be approximated for a practical implementation 
here. The first approximation is given by [95] 


D i:j = 0.002628 


mg) 


(55) 


where the units of Dij is err? Is, p is atm, T is K. 

The reduced mass and combined Lennard-Jones (12-6) parameters contained in the 
preceding equations are each computed according to the equations 


rriimj 
rrii + m,j 

(56) 


(57) 

II 

rbl 

<S>. 

Cb 

| 

(58) 
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5.2.4 Specific Heat, Enthalpy, Entropy, and Gibbs free energy 

The specific heat at constant pressure, enthalpy, and entropy were computed for each species 
in the gas mixture according to the polynomial equation used by McBride, et al. [119], given 
as 

(j 

= d\T 2 + cl 2 T 1 + a 3 + q 4 T + + clqT^ -j- a 7 T^ (59) 

h T T 2 T ^ T ^ 

—— = -aiT” 2 + a 2 T~ 1 In T + a 3 + a 4 — + a 5 — + a 6 — + a 7 — + a 8 T _1 (60) 

lly 1 ^ O 4 U 

g Jp—2 rp 2 rpS rp4 

— = -ai— a 2 T~ 1 + a 3 InT + a 4 T + a 5 — + a 6 — + a 7 — + a 9 (61) 

iiy ^ ^ o 4 

where the coefficients for each species are given in the same reference and is the some form 
used by the NASA Glenn program CEA. There are two sets of coefficients for each equation, 
one for the temperature ranging between 200 K and 1000 K and another for the temperature 
ranging between 1000 K and 6000 K. Computing the thermodynamic properties of the flow 
in such a matter is a simple yet accurate method to do so. 

The computation of Gibbs free energy, required when computing the equilibrium con- 
stant of a reaction, uses the thermodynamic relation 

g = h — Ts (62) 


where enthalpy and entropy are simply the values computed from Equations 60 and 61. 

5.3 Solution Methodology 
5.3.1 MacCormack’s Method 


MacCormack’s method [120] is a finite difference predictor-corrector scheme that is second- 
order accurate in both time and space. The method has two steps. In the first step the 
conserved variables in the temporal term are calculated at a temporary time step by the 
equation 

Upf = u * + At(^y (63) 

Note that the temporary time step is denoted by t + At. The time derivative in Equation 
63 is given by 




(64) 
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where the spatial derivative is given by a first-order forward difference at the initial time 


step. At this point the elements of the conserved vector U- +A< may be used to obtain the 
primitive variables for use in the corrector step. 

The corrector step obtains the conserved variables at the next time step according to 
the equation 

... 1 r — ^ /r)TT\ t+At 

(65) 


u: 


■t+At 


u: 


t+At 


+ u + At — 


MJi 


Again, the time derivative at the temporary time step is given by 


m\ t+At __ftm t+At _ n t+A t 
dt)i V dnj , 


( 66 ) 


where the spatial derivative is given by a first-order backward difference at the temporary 
time step. This process is repeated until a convergence criterion is satisfied. 

5.3.2 Numerical Damping 

Numerical damping is required in a shock-capturing scheme in order to dampen the Gibbs 
ringing phenomenon around areas of large gradient in the solution domain such as the shock 
wave. The term describing numerical damping at each grid point has the form 


where 


S<j> \4>i+i - 24>i + (j)i-i\ dU 

D <t> = 7T vJ ; iirK + a i)^T~ 

(1 - U^cpTVD + on 


(f>TVD = \4>i+l + 4>i\ + 1 4>i — 01-1 1 


(67) 


(68) 


— 4 > i + 1 + 20i + <t>i - 1 


(69) 


a i = \ 7 


Pi 
/ — 

Pi 


(70) 


The placeholder variable 0 represents any variable used in the damping expression, for 
example pressure, temperature, or species concentration. If multiple variables are used in 
the damping scheme simply add each to the PDE, which at spatial step i takes the form 

t / rnv T X t 


(W); + + Q * + (D 01 + D ^ + 


(71) 
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Figure 127: Comparison of the CHEMEQ algorithm and the CHEMEQ algorithm using 
a-QSS. [84] 

The damping term should be updated for each step of the MacCormack method. Note 
that when differencing the flux vector in the damping equation, the opposite type of dif- 
ferencing should be used. For example, if forward differencing is currently being used in 
the MacCormack scheme, backwards differencing should be used for the damping term and 
vice versa. 

5.3.3 Finite-Rate Chemistry 

The chemistry solution uses an algorithm described by Mott, et. al [84] known as a-QSS. 
The a-QSS algorithm is a methodology designed specifically for solving the stiff system of 
equations present in reacting systems. In his paper, Mott demonstrates that when the com- 
monplace CHEMEQ algorithm’s hybrid method 2 is replaced by the a-QSS method (called 
CHEMEQ2), more rapid convergence is achieved over standard CHEMEQ as shown in Fig- 
ure 127. The algorithm also benefits from very small overhead and storage requirements. 

The method begins with the decomposition of the species rate of change term into 

2 A hybrid method is one that chooses stiff or non-stiff integration depending on the time scale of the 
equation to be integrated. 


167 


production and loss terms by the equation 


dcj 

dt 


Gi - piCi 


(72) 


where Gi is the production rate of species i and ptCi is the loss rate. This form of the species 
production equation is chosen because for constant p and G, the equation has an exact 
solution that forms the basis of this method (and quasi-steady-state methods in general). 
In fact, for G and p which are slowly varying over a particular timestep, the exact solution 
represents an approximate, first-order version of a quasi-steady-state algorithm. 

The exact solution of Equation 72 is recast in an convenient algebraic form to yield 
the basis for a predictor-corrector formulation upon which the algorithm depends. The 
predictor and corrector equations are simply 


C„ = Co + Atchern ^7 (73) 

p chem l + aAt chemPo K > 

Cc = Co T At c hemZ : VT — (74) 

1 + aAt chem p 

where the subscript 0 indicates the initial value, p the predicted value, and c the corrected 
value. The predictor-corrector equation set is iterated upon (with the corrected value of 
one iteration being used as the predicted value of the next) until convergence is achieved. 
Weighted parameters appearing in Equations 73 and 74 are given by 


P = ^(PO+Pp) 


G = aG p + (1 — a)Go 


with 


L, 


Po = 


Pp = 


0 


Co + £r, 


Lr 


Cp T fy 


(75) 

(76) 

(77) 

(78) 


The inclusion of the machine epsilon terms in the denominators is a practical consideration 
in order to avoid divide-by-zero errors for species with no concentration. 

The a parameter, which arises due to the convenient algebraic form of Equation 73 and 
74, is actually a function of pAt and contains exponential terms. Because exponentials 
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are not well suited to numerical computation, a polynomial approximation to a is desired. 
Using the Pade approximation, 


180r 3 + 60r 2 + Hr + 1 
360r 3 + 60r 2 + 12r + 1 

180r 3 + 60r 2 + llr + 1 
360r 3 + 60r 2 + 12r + 1 


where 


f — (PO^tchem "P ^m) 


(79) 

(80) 


(81) 


r = (pAt chem + e m ) 1 


(82) 


If at any point during the iteration cycle c c — c p < e m c c for all species, the chemistry is 
converged over the chemical time step for the cell. 

While it is tempting to simply run the ce-QSS predictor-corrector algorithm over the 
entire fluid dynamic time step, stability considerations dictate that the fluid dynamic time 
step actually be partitioned into much smaller chemical time steps with chemical integra- 
tions occurring over each partition. The initial chemical timestep is chosen to be the same 
as for the previous global iteration (with the initial timestep chosen to be some arbitrarily 
small value). If the predictor-corrector loop fails over a particular chemical timestep, the 
iteration is restarted by multiplying the chemical timestep by 1/3. If the predictor-corrector 
step is successful, the time step is allowed to grow according to the equation 


At c hern,new chem.old ( T 0.005 

' a 


(83) 


where 


v Vc 

with v = 0.001 and a chosen to be the largest among the different species. 

5-4 Validation 


(84) 


The ability of QUIVR to properly predict the physics of high temperature, reacting flow 
along the stagnation streamline will now be summarized. Validation can be achieved by 


169 


comparing the results produced by QUIVR against experimental, theoretical, and numerical 
results present in the literature. The validation questions can be summarized as follows: 


Can QUIVR... 

• ...successfully reproduce the stagnation point heat transfer experienced by stagnated 
flows? 

• ...reproduce the space- varying thermodynamic properties in the shock layer? 

• ...match the reaction front oscillation frequencies and regimes of stability observed in 
the literature? 

5.4.1 Stagnation Point Heat Transfer 


Comparing the wall heating computed by QUIVR to the theoretical relationship of Fay and 
Riddell [121] is a useful validation study to ensure the code is properly computing heating 
and boundary layer properties. The Fay and Riddell relation for wall heat flux without 
chemical dissociation is given by the equation 


q w = 0.76 Pr U ' b (p e /u) a4 (p^ J ai \/^(ho - M 


ds 


l + (Le° -52.!)^ 
ho . 


(85) 


where the overline above certain quantities indicate that the quantity is dimensioned. Fur- 
thermore, the subscript e denotes boundary layer edge properties and w denotes wall prop- 
erties. Because the velocity vector is directed toward the body, the boundary layer edge is 
defined here to be the distance from the wall where the stagnation enthalpy is within 1% 
of the freestream value. It should also be noted that the numerical correlation of Equation 
85 is only valid for a range of parameters for Pr = 0.71, those being 


• Le = 1.0, 1.4, or 2.0 

• 300 K < T w < 3000 K 

• 1770 m/s < Voo < 6950 m/s 

• 7620 m < h < 35,570 m 
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Because the QUIVR solution is in terms of dimensionless quantities, it is convenient to 
recast the Fay and Riddell relation in a nondimensional form. This is accomplished simply 
by substituting each dimensional term in the equation by the nondimensional version times 
the nondimensionalizing factor. This results in the following equation, with the dissociative 
enthalpy term removed because this validation study is not considering dissociation, 

q w = °^— Pr~ 0 - 6 (p e p e ) 0A (p w p w ) 0 - 1 J d ^(h 0 - h w ) (86) 

y/Reoo V as 

Note that the du e /ds term references the edge velocity gradient in the direction paral- 
lel to the surface. In most applications it is typical to use the approximation du e /ds ~ 
l/ r n\/2{Pw ~ Poo)/ Pw i which is the modified Newtonian approximation of the edge veloc- 
ity gradient. However, because the quasi- ID formulation actually solves this term directly 
(recall that u\ = du/ds ), it is best in this case to simply utilize the solved value of u. 

The test case used to validate the heat flux calculation was selected to he within the 
valid range of the Fay and Riddell relation. A 63.5 mm diameter projectile traveling at 1770 
m/s in ambient air at 35 km conditions and with a wall temperature of 3000 K was chosen 
as the test case. At 35 km, Too = 237 K, = 558 Pa, and Moo = 5.72 (from the velocity 
chosen), and Re oo = 29,470. The gas was chosen to be calorically perfect throughout the 
flow. The 401 point computational grid with n max = 0.18 was clustered near the wall in 
order to get high boundary layer resolution as heat flux is highly sensitive to spacing near 
the wall. 

The converged heat flux as a function of distance from the projectile wall is shown in 
Figure 128. Note that in the freestream and behind the shock wave the heat flux remained 
at zero until the boundary layer is encountered. There is a nonphysical spike in the heat 
flux going though the shock wave which is a consequence of the shock capturing nature of 
the solution. The width of the spike, corresponding to the thickness of the shock wave, is of 
little consequence to the solution sought here where high resolution near the wall is preferred 
to a thinly resolved shock wave. Table 21 shows a comparison of the QUIVR solution of 
wall heat flux and the value obtained from Equation 86 using edge and wall values obtained 
from the QUIVR solution. There was observed to be only a 0.34% difference in wall heating 
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Figure 128: The solution of nondimensional heat flux as a function of distance from the 
wall is shown for a 63.5 mm projectile with T w = 3000 K and traveling at = 5.72 at 
an altitude of 35 km. 

Table 21: Comparison of Fay and Riddell with the QUIVR solution of heat flux for the 
validation case. 



(Pl)e 

{pl)e 

(pi)w 

(Ml )w 

Ol)e 

(^l)o 

(^1 )w 

(Q2)w 

QUIVR Solution 

5.584 

4.230 

3.363 

5.941 

0.718 

0.576 

0.966 

0.00855 

Fay and Riddell 

- 

- 

- 

- 

- 

- 

- 

0.00858 

% Difference 








0.34 % 


values, thereby lending confidence to the handling of heat flux in the QUIVR code. 

5.4.2 Thermodynamic Properties in the Shock Layer 

Tracking the variation of properties in the shock layer of a sphere in hypersonic flow is 
another way to validate the behavior of the QUIVR code. Freestream conditions were ref- 
erenced from Kumar and Graves [122], which provides a data set against which to compare 
the results of QUIVR. While the reference data was computed for a spherically-blunted 
cone, along the stagnation streamline the solution was comparable to the spherical case. 
The simulated freestream conditions were = 10.33, p Q 0 = 100.77 Pa, and T , = 46.26 K. 
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Figure 129: Comparison of QUIVR p/p s hock result in shock layer compared to reference 
data. 

These conditions approximately correspond to high-altitude, hypersonic flight. The simula- 
tion used a thermally perfect mixture of oxygen and nitrogen at atmospheric proportions. 
The wall temperature was set isothermally as 330.6 K. The nose radius was set to 0.03175 
m. The grid consisted of 401 grid points equally spaced from the wall to n max = 0.15. 

The results of the analysis in Figures 129 and 130 show p/p s hock and T /T s h oc k as a 
function of distance from the wall normalized by shock standoff distance. It was observed 
that the numerical solution can accurately reproduce the reference density and temperature 
profiles from the shock to the wall. It was also observed that the boundary layer thickness 
was well captured by the simulation. 

5.4.3 Qualitative Comparison to Experiment and Numerical Simulation 

The experiments of Ruegg and Dorsey [96] were simulated using the quasi-ID code in order 
to qualitatively assess the suitability of the proposed analysis method. Ruegg and Dorsey 
conducted early experiments on PMSIC by firing spherical, 20 mm diameter projectiles into 
a stoichiometric hydrogen-air mixture in the vicinity of Mach 5 at pressures ranging from 0.1 
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n /$s 


Figure 130: Comparison of QUIVR T/T s } loc h result in shock layer compared to reference 
data. 


to 0.5 atm. Over the range of experimental conditions investigated, three distinct regimes of 
behavior of the reaction/shock front system were observed: stable; regular, unsteady; and 
large-disturbance, unsteady. A qualitative comparison between the quasi- ID simulation 
results and experimental shadowgraphs at the conditions in Table 22 will now be made. 
The simulated results consist of an evolution of water mass fraction and pressure on the 
stagnation streamline as a function of non-dimensionalized time and distance (n) from the 
surface of the projectile. 


Table 22: Selected Ruegg and Dorsey experimental conditions for a 20 mm diameter, 
spherical projectile fired into premixed, stoichiometric hydrogen/air. 


Case 

Moo 

Poo (atm) 

Too (K) 

Stability Regime 

1 

5 

0.10 

300 

Stable 

2 

4.9 

0.25 

300 

Regular, Unsteady 

3 

4.3 

0.5 

300 

Large Disturbance 


Figure 131 shows quasi-ID simulation and experimental results for the stable regime 
of oscillation. The simulated and experimental results indicate reaction and shock front 
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locations that were invariant with time. This time independence was observed experimen- 


tally by smooth shock and reaction fronts. The reaction and shock fronts were separated 


by a relatively large induction distance with no direct interaction between the two beyond 


the shock wave raising the temperature of the hydrogen/air mixture above its autoignition 


limit. There was qualitative agreement between simulation and experiment for this case. 
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Figure 131: (a) Contours of water mass fraction on the stagnation streamline overlaid 
with contour lines of pressure as a function of non-dimensional time for Ruegg and Dorsey 
experiment Case 1. (b) An annotated experimental shadowgraph captured by Ruegg and 
Dorsey for the same conditions is shown for comparison. 


The regular instability regime is reproduced in Figure 132. The experimental shadow- 
graph shows a corrugated pattern in the wake of the sphere caused by an oscillation in 
density due to the periodic creation of new reaction fronts in the stagnation region of the 
flow. Between the corrugated, burned region of gas and the shock wave lies a region of un- 
burned gas. This is due to the induction time between the bow shock and the reaction front. 
Lines of density variation in the unburned region are due to periodic contact discontinuities 
traveling from the bow shock to the beginning of each new reaction zone. The contact 
discontinuities are not visible in the simulated results in this figure because contour lines of 
pressure are shown (pressure being equal on both sides of a contact discontinuity). The con- 
tact discontinuities in the computations will be visible in later figures when the computed 
results are compared to established ID mechanisms. There is qualitative agreement between 
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the experiment and the present computations because the physical mechanism observed in 
the simulated results match the flow features witnessed in the experimental shadowgraph. 


Reaction Bow 

Front Shock 


0 0.05 0.1 0.15 

n 

(a) 

Figure 132: (a) Contours of water mass fraction on the stagnation streamline overlaid 
with contour lines of pressure as a function of non-dimensional time for Ruegg and Dorsey 
experiment Case 2. The drawn contact discontinuity is approximated based on separate 
density contours, (b) An annotated experimental shadowgraph captured by Ruegg and 
Dorsey for the same conditions is shown for comparison. 




(b) 


The large-disturbance instability regime is shown in Figure 133. Experimental results 
of the flow at the Case 3 condition show a large-amplitude heaving of both the bow shock 
and the reaction front. Both the bow shock and reaction zone oscillate with a much lower 
frequency than that observed in the regular regime and apparently non-periodic features 
are visible in the shadowgraph. Although both the experiment and the present numerical 
results show large amplitude behavior, it is difficult to isolate other features to determine 
how well the quasi- ID simulation is performing. For further comparison, the present results 
are compared to a higher-fidelity, 2D simulation. Figure 134 shows a zoomed-in view of a 
section of the present quasi- ID simulation compared to a segment of a simulation at the same 
conditions conducted by the author using VULCAN. These images show contours of density 
as a function of time and distance along the stagnation streamline. The present results show 
a striking agreement with the high-fidelity numerical result and clearly show coinciding 
features of reaction shocks, retonations, contact discontinuities, and a superdetonation as 
the bow shock is penetrated by the reaction front. These features will be discussed in further 
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detail later, but for now it appears that the quasi-ID results are capable of predicting the 
large-disturbance instability regime. 
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Figure 133: (a) Contours of water mass fraction on the stagnation streamline overlaid 
with contour lines of pressure as a function of non-dimensional time for Ruegg and Dorsey 
experiment Case 3. (b) An annotated experimental shadowgraph captured by Ruegg and 
Dorsey for the same conditions is shown for comparison. 


5.4.4 Quantitative Reaction Front Oscillation Frequency 

The current analysis approach has been applied to simulate the ballistic experiments con- 
ducted by H.F. Lehr [10] in which 15 mm diameter projectiles were fired into a stationary, 
premixed, stoichiometric hydrogen/air mixture at near-hypersonic speeds. The ambient 
conditions of the fuel-air mixture had a static temperature and pressure of 292 K and 0.42 
atm, respectively. Before the results from the current study are discussed, some quick ob- 
servations are made about Lehrs experiments and the numerical simulations by Sussman 

[98]- 

Lehr investigated the behavior of the reacting system over a wide range of velocities, 
1685 to 2058 m/s. For each case, the frequency of oscillation of the reaction front was 
determined using shadowgraph images and considering the distance between successive 
striations in the wake of the sphere. Lehr found that as projectile velocity increased the 
frequency of oscillation increased monotonically. At high projectile velocities an overdriven 
detonation occurred which caused an apparently stable reaction front. 
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Figure 134: A detail view of density contours on the stagnation streamline as a function of 
time for Case 3 for the (a) QUIVR quasi- ID simulation and (b) high-fidelity, 2D simulation 
produced using VULCAN. Note that the computed results are non-dimensional in time 
while the high-fidelity plot is in seconds. 

The Chapman- Jouguet detonation velocity of a stoichiometric hydrogen-air mixture 
with initial conditions as that of the Lehr experiments was calculated to be 1957 m/s using 
the Chemical Equilibrium with Applications (CEA) program developed at NASA Glenn 
by McBride and Gordon [123, 114]. This velocity is sensitive to species present — in this 
case only the species H, O, OH, H 2 O, HO 2 , H 2 , O 2 , and N 2 were allowed as product 
species. The experimental value reported by Lehr was 2030 m/s, corresponding to the 
velocity past which oscillations in the reaction front were no longer visible. It should be 
noted, however, that other analyses have established that oscillations might occur above the 
Chapman-Jouguet detonation velocity. For example, the experiments of Alpert and Toong 
[12] found oscillations in the overdriven range of flight and modern numerical analyses have 
also substantiated this observation [110, 111, 112]. 

Sussman [98] simulated the axisymmetric flow around a 15 mm diameter sphere at 
the same ambient conditions and freestream velocities studied by Lehr using an upwind, 
TVD scheme. The oncoming flow is composed of a stoichiometric hydrogen-air mixture. In 
general, the numerical results of Sussman matched well with the experimental results with 
the exception of cases near the high velocity limit. A large discrepancy was observed at 
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2058 m/s, which was postulated to be due to the method of interpreting the experimental 
frequency at that point. However, the experimental oscillation frequency at this point 
was successfully replicated in numerical work done around the same time by Ahuja [110, 
111, 112]. At a projectile velocity of 2119 m/s, Lehr showed no oscillation although the 
numerical work by Ahuja showed oscillations occurring all the way to 2605 m/s; however, 
at this point they were of very small amplitude and located away from the stagnation 
streamline. Sussman was able to simulate oscillations for the 2119 m/s case but not for the 
2605 m/s case. 

Figure 135 and Table 23 show results from the quasi-ID simulations over the velocity 
range studied by Lehr and Sussman. Note that due to the reduced computing requirements 
of the quasi-ID simulation, a large number of cases are computed over the velocity range 
of the experiment versus previous numerical studies. In doing so, a better idea as to the 
general shape and limits of the parameter space is obtained. It is observed that, in general, 
the present results are in reasonable agreement with experimental results and higher-order 
simulations. Overall the quasi-ID simulations over-predict oscillation frequencies by ap- 
proximately 20%-30%. In all cases the overall behavior of the flow is predicted (e.g. regime 
of instability) well. 


Table 23: Oscillation frequencies from experiment and numerical investigations. 


Velocity (m/s) 

Lehr (kHz) 

Sussman (kHz) 

Present Results (kHz) 

1685 

148 

135 

185 

1804 

425 

420/74 

- / 98 

1931 

712 

711/103 

881/- 

2029 

1040 

1038 

1354 

2058 

1960 

1180 

1423 

2119 

n/a 

1500 

1851 

2605 

n/a 

n/a 

n/a 


One of the benefits of having the ability to more densely populate the projectile velocity 
parameter space is that transition points from one regime of stability to the next become 
more clear and the general morphology of the space begins to focus. In addition, transition 
points may be more finely interrogated by the quasi-ID simulation in order to observe 
the changing behavior of the flow in order to arrive at mechanisms for transition. With a 


179 


2000 


o Lehr (Experiment) 
a Sussman (Axi-2D Simulation) 
■ Axdahl (Quasi-ID Simulation) 


1500 


N 

X 

> 

u 

c 

cu 

=5 

cr 

cu 


1000 


500 


Aoo 


1600 


1700 1800 1900 2000 

Projectile velocity (m/s) 


2100 


2200 


Figure 135: Reaction front oscillation frequency as a function of projectile velocity for 
experimental and numerical investigations. Note that frequencies are grouped due to FFT 
windowing. 
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more densely populated parameter space, trends not observed in previous work also appear. 
For example, the space takes on something akin to a terraced appearance where there are 
distinct regions of stability with rapid transitions in some locations from one dominant 
frequency to the next (e.g. around 2000 m/s and 2080 m/s). 

The basic transition mechanism at work in the observed results is an out-of-phase be- 
havior of pressure waves generated by the new reaction fronts. For projectile velocities 
below 1650 m/s, the reaction front is observed to be smooth. If new reaction fronts are 
being created, they are imperceptibly low in amplitude and strength at the studied grid 
resolution. At the beginning of the regular regime zone, there are moderate amplitude 
reaction front features. Increasing the projectile velocity has the effect of creating addi- 
tional reaction zones along the same reaction peninsula. This is due to the pressure waves 
generated by new reaction zones being significantly out-of-phase with the reflected pressure 
waves from the projectile surface. As velocity increases, new reaction fronts are created 
by out-of-phase waves propagating faster than expansion waves can relax the shock front. 
This causes a secondary, low-frequency oscillation to appear mainly between 1730 m/s and 
1960 m/s upon which the higher-frequency regular instability is superimposed. The sec- 
ondary oscillation dies out around the detonation velocity of the mixture (1957 m/s), but 
this may be coincidental. As the projectile velocity increases past 1960 m/s, where the 
classic, in-phase regular regime is re-established, a terraced trend of oscillation frequency 
versus projectile velocity occurs. The cause of this trend is again due to narrow velocity 
ranges where out-of-phase pressure wave behavior creates new reaction fronts in addition 
to ones already being created periodically due to the classic regular regime mechanism. 
The reaction front at the transition points shows mixed regular regime behavior with non- 
cohesive, irregular reaction fronts created by the out-of-phase pressure waves. In general, 
as the projectile velocity increases, the amplitude of the reaction front is found to decrease, 
eventually smoothing out into apparently stable behavior as it begins to couple with the 
shock front. 
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5.4.5 Comparison to ID Mechanism Models 


The simulated experiments of Ruegg and Dorsey were reexamined in order to compare the 
present quasi-ID numerical results to published ID instability mechanism models in the 
literature. First, the quasi-ID simulation of the regular, unsteady regime was compared 
with the published mechanism of Matsuo [124] in Figure 136. Note that all of the major 
features of the mechanism proposed by Matsuo are present in the simulated result. A newly 
created reaction front produces pressure waves that travel downstream toward the projectile 
surface and upstream toward the bow shock. The upstream traveling pressure wave meets 
the bow shock roughly in phase with the pressure wave reflected from the projectile surface. 
The collision of the pressure waves and the shock front causes the shock to move upstream 
with a contact discontinuity reflected downstream. The contact discontinuity is visible on 
the density plot and is visible as a line traced from each cusp of the bow shock to the 
newly created reaction front. The hot gas upstream of the contact discontinuity has a lower 
induction time than the downstream gas, which creates a new reaction zone that eventually 
extinguishes the old one. This process repeats indefinitely for a projectile with constant 
velocity. The simulated result differs slightly from the published mechanism in that the 
simulated contact discontinuity is connected to the tip of the reaction zone peninsula versus 
the ID mechanism, which traces it to the point where the old reaction front is merged with 
the new one. In addition, the presently simulated result indicated that the reflected pressure 
waves skip a generation of new reaction zone before interacting with the bow shock. The 
published mechanism shows the reflected pressure wave interacting with the bow shock along 
with the upstream-traveling pressure wave of the next reaction zone generation. Finally, the 
present computations showed a rarefaction wave reflected at the point where the pressure 
wave and the bow shock meet, which is not depicted in the published model. 

The quasi- ID simulation of the large-disturbance, unsteady regime, shown with a density 
contour plot, will now be compared with the corresponding mechanism of Matsuo and Fujii 
[97] in Figure 137. The large-disturbance mechanism of Matsuo and Fujii requires an energy 
release large enough to create a deflagration to detonation transition (DDT), noted in the 
figure as an explosion in an explosion (this is discussed further in the detonation literature 
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Figure 136: Contour lines of (a) density on the stagnation streamline as a function of time 
for Ruegg and Dorsey experiment Case 2 and (b) ID mechanism for the regular, unsteady 
regime by Matsuo (reproduced with permission). 


[106, 107]. This causes a penetration of the self-propelled reaction front through the bow 
shock, reflecting a rarefaction wave and contact discontinuity downstream at the point of 
intersection. The self-propelled reaction front, which is coupled with the bow shock for 
a short time, gradually separates from the shock wave and relaxes to its original position 
before resuming the cycle. Comparing the present results to the ID mechanism shows that 
all of the essential features of the mechanism are reproduced in the quasi- ID simulation. 
Figure 10 shows a detailed view of the detonation region of the quasi-ID simulation to 
examine what is happening closer to the shock wave penetration. Here it is clearly visible 
that the intersection of the reaction shock wave from the first explosion with the bow shock 
is responsible for generating the rarefaction wave and the contact discontinuity that causes 
the second explosion. The second pressure wave seen in the simulated results is due to 
the generation of the second explosion, explaining why it doesnt appear in the general, 
single-explosion mechanism. After the second explosion, the pressure wave and reaction 
front traveling upstream merge with the bow shock and a self-sustained superdetonation 
is observed for a short amount of time. Overall the quasi-ID simulation shows excellent 
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agreement with the ID mechanism found in the literature. 


a> 

E 


Reaction Bow 
Front Shock 



0.05 


0.15 

n 


i . 

0.25 



(a) 


(b) 


Figure 137: Contour lines of (a) density on the stagnation streamline as a function of time 
for Ruegg and Dorsey experiment Case 3 and (b) ID mechanism for the regular, unsteady 
regime by Matsuo and Fujii (reproduced with permission). 


5. 5 Effect of Equivalence Ratio on Stability 

It has been shown that QUIVR was able to reproduce the important physics occurring in 
blunt-body PMSIC at a much lower computational cost versus axisymmetric simulation. 
One benefit of this capability is the ability of the computational framework to conduct 
parametric simulations of the bring of a blunted projectile at different conditions. Some 
attempt has been made to do this in the literature, but knowledge gained is limited due 
to the computational time required for multidimensional simulations to converge as well as 
the post-processing of numerical data. 

For the purposes of this study, QUIVR was integrated into DAKOTA 3 , an automated 
code execution and analysis framework, in order to allow the automated execution of the 
program over the design space of interest. Here the design space is bounded by projectile 
velocities ranging from 1650 to 2500 m/s and ambient pressures ranging from 0.05 to 1.0 

3 Code obtainable from http://dakota.sandia.gov/ 
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Figure 138: A detail view of the (a) stagnation streamline density contours from the present 
simulation (Case 3) compared to a (b) ID sub-mechanism of the large-disturbance regime 
as proposed by Matsuo and Fujii [97] (reproduced with permission). 


atm. The design space is partitioned into 256 equally-spaced cases and the entire set of cases 
were run repeatedly for equivalence ratios of 0.5, 1.0, and 1.5. The ambient temperature 
was held constant at 292 K and the projectile diameter was set to 15 mm, equivalent to the 
projectile used in the Lehr study. In order to facilitate the post-processing of the data, the 
wall water mass fraction is tracked as a function of time and then analyzed with a discrete 
Fourier transform automatically at the conclusion of each simulation. The chemical model 
was chosen to be the Jachimowski reduced 8-species, 13-reaction model used in Chapter 3. 

Figure 139 shows filled contours of reaction front frequency as a function of ambient 
pressure and projectile velocity for the lean, stoichiometric, and rich cases. These plots are 
unique contributions to the literature and are a direct result of the low computational cost 
to run the QUIVR code in a reacting mode. Each space is characterized by an unstable 
region with reaction front frequencies in excess of zero kilo-Hertz bracketed on either side by 
stable regimes. To the left of the unstable region are stable cases where the reaction front 
is detached from the shock front and to the right is a region where the reaction and shock 
fronts have merged into a stable detonation. The region of regular instability shows areas 
where the reaction front frequency decreases and then increases again. Such areas are due 
to the envelope regime of stability where the larger amplitude, lower frequency disturbance 
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is detected by the FFT. 

As the equivalence ratio is increased the regime of instability shifts to higher velocities. 
This is a Mach number effect in that increasing the equivalence ratio of a hydrogen-air 
mixture decreases the mixture molecular weight. Because for a particular case in the pa- 
rameter space temperature, pressure, and velocity are being held constant, the freestream 
Mach number of the flow increases with decreasing gas constant. Increasing the flow Mach 
number caused a higher post-shock temperature that in turn decreased the induction time, 
shifting the behavior of the reaction front to a higher effective velocity relative to the base- 
line stoichiometric case. Note that these trends would presumably be reversed — increasing 
equivalence ratio would shift the trends to the left — for a heavier-than-air fuel as increasing 
equivalence ratio in that case would increase the molecular weight of the mixture. 

Figures 139(b) and 139(c) both indicate the Chapman-Jouget detonation velocity through 
the parameter space, solved using CEA. The CJ line appears as a vertical line because it is 
only a weak function of pressure for the conditions considered. It was shown that instabilities 
may occur above the Chapman-Jouget velocity, which has been observed in experimental 
and numerical literature and was confirmed here. 

5.6 Chapter Summary 

The primary concern of this chapter was the development of a numerical methodology for the 
rapid solution of the premixed, shock-induced combustion ahead of blunted projectiles. Such 
flows are characterized by steady or unsteady physics depending on the conditions of the 
oncoming fuel-air stream. The current paradigm for studying such flows is an axisymmetric 
simulation in the vicinity of the stagnation region, but such simulations are still taxing due 
to the requirement for fine grids and the solution of chemical reactions using global time 
stepping. The currently proposed method, a quasi-ID form of the axisymmetric Navier 
Stokes equations, considerably reduces the computational requirements to solve such flows 
and makes such solutions practical on a personal computer or high performance workstation. 

It was shown that the quasi- ID framework, implemented through the QUIVR code, 
was able to successfully reproduce the important physics characteristic of shock-reaction 
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Figure 139: Filled contours of reaction front frequency as a function of projectile velocity 
and ambient pressure for (a) lean, (b) stoichiometric, and (c) rich equivalence ratios. The 
black line on the latter two plots indicate the Chapman-Jouget detonation velocity. 
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front instabilities at selected conditions. This supports past experimental and numerical 
observation that the important physics of such flows originate in the stagnation region 
of blunted bodies. For the regular regime of instability, the important mechanism is the 
interaction of pressure waves emanating from newly created reaction zones with the bow 
shock, reflecting contact discontinuities that create further new reaction zones. For the large 
disturbance regime, the important mechanism is the creation of a deflagration-to-detonation 
transition that causes a rapid upstream advance of the reaction front that penetrates the 
bow shock. These mechanisms were captured by the quasi- ID formulation. 

Past numerical studies have attempted to parametrically study PMSIC due to blunted 
projectiles but were limited by the computational demands of an axisymmetric simulation 
and the post-processing required for each case. Because of the significantly reduced demands 
of the current framework, such studies are able to be rapidly conducted. One example study 
is the effect of fuel-air equivalence ratio over given ranges of projectile velocity and ambient 
pressure. The current study shows results on these effects that have not existed in the 
literature to this point. 
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CHAPTER VI 


SUMMARY AND AVENUES OF FUTURE WORK 

6.1 Summary of Contributions 

The application of premixed, shock-induced combustion to airbreathing engines has been 
identified as an enhancing technology to access the hypervelocity regime of flight. By mixing 
fuel with air on the forebody of a shcramjet vehicle, the internal flowpath need only provide 
for combustion of the fuel-air mixture. This allows for significantly shortened vehicles over 
those that require injection, mixing, and combustion to occur internally. 

Previous studies have identified three main issues with the use of PMSIC for high speed, 
airbreathing propulsion. The first is the ability to inject and mix fuel with a high Mach 
number airstream while sustaining the overall performance of the flow. The second is 
the prevention or mitigation of autoignition that may occur on the forebody as the fuel- 
air mixture accesses hot areas of the forebody flow, such as the boundary layer. The 
third is the possible presence of instability in the shock-reaction front at the entrance 
of the combustor as the literature shows that shock-reaction front interactions may give 
rise to either complete stable behavior or highly unstable behavior. It was the goal of 
this dissertation to systematically address each of these issues via studies using numerical 
simulation. 

The study of mixing and injection at conditions relevant to Mach 12 flight yielded new 
analysis capabilities and contributions to knowledge. Hydrogen injector geometries within 
the flush-wall, ramp, and strut injector classes were simulated and compared at the same 
freestream conditions. A new method compared each injector on the basis of the losses 
incurred for any given level of mixing achieved. This allowed the technically defensible 
identification of injector concepts that will best meet any given mission needs. It was shown 
that the strut injector was able to achieve high levels of mixing with a moderate amount 
of stream thrust potential loss. The ramp injector was only able to achieve modest mixing 
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gains, but was able to preserve the highest level of stream thrust potential at cowl closure. 
The flush- wall injector was able to achieve similar levels of mixing as the strut but at inferior 
levels of stream thrust potential. The decision of what injector to use for any given flowpath 
is dependant on the respective system requirements on losses allowed and required level of 
mixing, although it appears that, from a performance standpoint, flush-wall injection is at 
a disadvantage for the geometries and conditions studied. The usefulness of the flush- wall 
injection strategy for forebody injection was further weakened after conducting reacting 
simulations that showed extensive burning in the shock layer ahead of the fuel plume. 

The mitigation of autoignition at the wall was addressed through the use of wall cooling. 
After studying the parametric effect of cooling geometry, angle of injection, species, coolant 
temperature, and coolant mass flow rate, it was learned that hydrogen and helium coolants 
were both well suited for mitigating autoignition when injected tangentially to the wall 
through film cooling slots. The use of hydrogen as a film coolant is beneficial to a flight 
system due to the coolant being the same species as the fuel, although the levels of coolant 
mass flow rate required for substantial autoignition mitigation may pose an issue for vehicle 
propulsion efficiency. This drawback may be mitigated by considering a fueling split between 
the film cooler and the injector, however, to target a unity global equivalence ratio. It was 
also learned that nitrogen (and, similarly, other high molecular weight, inert species) are 
inappropriate for wall cooling when autoignition mitigation is required due to its poor 
specific heat by mass. In addition, it was found that transpiration cooling was unable to 
mitigate wall autoignition to any useful degree for the studied conditions and configurations. 

Finally, the issue of reaction front/shock wave instability was studied through the rep- 
resentative problem of premixed, shock-induced combustion ahead of a sphere traveling 
at hypersonic speeds. It was found that the axisymmetric Navier Stokes equations were 
able to be reduced to a quasi-ID form valid in the stagnation region of the sphere. By 
verifying the solution of these reduced equations against prior experimental and numerical 
results, it was found that the important unstable physics found in axisymmetric simulations 
were successfully reproduced. This confirmed that the important physics in reaction- front, 
boundary layer instability originate in the stagnation region of the sphere. The quasi- ID 
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equations were integrated into a tool and included in an automated analysis framework 
to enable rapid parameter studies. By doing such studies, it was possible to observe the 
regimes of instability in a multidimensional parameter space. This improved knowledge 
of the dependence of instability on velocity and pressure. For example, it was observed 
that there exists a single maximum for reaction front frequency and that the instability 
regime is encapsulated rather than extending unbounded in any given direction of pressure 
of velocity. 

The outcomes of the research presented in this dissertation both contribute to state- 
of-the-art knowledge of the physics of shock-induced combustion ramjets as well as define 
avenues for future study and development. Key results from this research indicate that 
forebody injection should be limited to ramp or strut injection depending on mission re- 
quirements constraining loss and required mixing. Autoignition should be approached using 
tangential film cooling using a light gas such as hydrogen or helium. Furthermore, instabil- 
ities arising at conditions relevant to the inlet of such a vehicle may be predicted using the 
proposed quasi-ID simulation framework. 

6.2 Avenues of Future Numerical Work 
6.2.1 Forebody Injector Design and Optimization 

The work presented in this dissertation was primary concerned with comparing different 
classes of injectors in the same type of forebody flow representative of hypervelocity flight 
and discussing the feasibility of each type of injector. The strut and ramp injectors both 
proved to be feasible for forebody fuel injection, although admittedly the designs used are 
point designs without regard to optimization. Parametrization and design space exploration 
of these geometries represent a future course of research with the current work forming a 
starting point. 

Ideally numerical optimization would take place on an actual forebody and would be 
automated to allow for more rapid design space exploration. Systematic exploration of a 
design space within an automated framework brings significant benefits in time, cost, and 
knowledge. The current design paradigm requires an iterative process that can be time 
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intensive and yield only point designs at each iteration. Design problems are compounded 
by the time issue in requiring a human in the loop at each discipline component of the design 
process. A paradigm shift in the current design process is obtained by integrating different 
codes in an automated framework. Such a shift yields significantly increased knowledge 
over reasonable variable ranges and allows the designer to achieve either optimized or robust 
designs through effective sampling of the design space and the application of state-of-the-art 
optimization routines. 

DAKOTA would be an ideal package to design this framework around. DAKOTA can 
accomplish parameter studies, optimization, or uncertainty analysis either by operating di- 
rectly on the design space or on surrogate models 1 constructed using design of experiments, 
an efficient sampling method also available in the toolkit. A further benefit of DAKOTA 
has the ability to change the type of study being conducted simply by changing the input 
file to the toolkit without further changing the overall framework. 

Figure 140 illustrates a notional N 2 diagram for a CFD framework that integrates a 
grid generator, numerical simulator, and post-processing package. DAKOTA is placed in 
the framework in a position where it can direct the code execution and flow of information 
between nodes in the diagram through the use of shell scripting. The generality of the nodes 
in the diagram are to emphasize that any code that has batch processing capability may be 
used in the framework. Examples of each type of program are indicated below the diagram. 

A DAKOTA-enabled CFD framework has the ability to systematically and automatically 
sample the design space using high-fidelity codes in order to construct fast-running surrogate 
models. A surrogate model is a polynomial equation that closely approximates the behavior 
of more complex, higher fidelity analyses over a defined design space. A typical surrogate 
model includes an intercept term, linear terms, interaction terms, and pure quadratic terms, 
as shown in Equation 87. 


n n — 1 n n 

R = Po + ^2 PiXi + 22 ^2 PijXiXj + 22 Pnrf 

i= 1 z=l j=i+ 1 z=l 


(87) 


1 Surrogate models are polynomial equations that closely approximate the behavior of complex, time- 

intensive codes over a given parameter space. 
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Figure 140: Notional N 2 diagram for the analysis framework along with example codes to 
be used in each node. 


After sampling the design space using parameter combinations defined from design of 
experiments in DAKOTA and constructing a surrogate model, the analyst may use the 
resulting model to further understand the design space through probabilistic design efforts. 

Design optimization is either carried out directly on a set of codes or on a surrogate 
model closely approximating the behavior of the design space. Optimization strategies that 
may be carried out in a DAKOTA-enabled framework include 


Gradient-Based Optimization Gradient-based optimization strategies may be 
carried out directly on a set of codes or on a surrogate model. In order to operate 
directly on a set if codes, the numerical simulator needs to be augmented with the 
ability to compute sensitivity derivatives in the neighborhood of a point design. In 
order to operate on a surrogate model, the design space needs to be representable by 
a polynomial equation (which may or not be the case depending on the problem being 
studied). Examples of gradient-based optimizers available in DAKOTA are Fletcher- 
Reeves conjugate gradient (CG), the method of feasible directions (MFD), sequential 
linear or quadratic programming (SLP or SQP), and the Broyden-Fletcher-Goldfarb- 
Shanno (BFGS) quasi-Newton technique. 
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Nongradient-Based Optimization For cases where the design space is not repre- 
sentable by a surrogate model and no derivatives are obtainable, DAKOTA may draw 
upon pattern search, Monte Carlo, and Solis- Wets algorithms. 

Stochastic Optimization DAKOTA has a highly-developed capability to use evo- 
lutionary (genetic) algorithms when operating on a set of codes or surrogate model. 
Single- and multi-objective genetic algorithms (SOGA and MOGA) are available for 
the analyst’s use, with the latter representing a powerful technique as most engineer- 
ing problems have multiple objectives. MOGA techniques do not require the use of 
weighted objective functions, instead seeking out points along the entire Pareto front 
of the design space. 

Design problems that may be encountered by the DAKOTA-enabled CFD framework 
can be addressed in many different ways without changing the structure of the framework 
itself. 


Large design spaces Design spaces that have a large number of parameters may 
be mitigated by employing a screening test before beginning higher fidelity, time- 
intensive analyses. A screening test is accomplished by sampling the design space using 
relatively few points (afforded by a two-level design of experiments) and investigating 
the results using an ANOVA. The ANOVA statistical technique will identify which 
variables produce the greatest impact on the variability of the response using the 
Pareto principle. Only the new, smaller subset of variables is used to define the 
design space with all others set at a default value. 

No derivative information If time or technical constraints preclude the use of 
gradient-based optimization techniques, nongradient-based techniques may instead be 
used in order to analyze the design space. For example, the use of an evolutionary 
algorithm will allow the identification of a global optimum without the use of any 
derivative information at any given point. 

Time-intensive codes If optimizing directly on a set of codes proves too time 
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intensive a task, the behavior of the design space can instead be defined using surrogate 
models. Surrogate models are constructed by systematically sampling the design space 
using design of experiments. Due to their polynomial nature, once a surrogate model is 
created the design space can be explored in real-time. Assuming it has been correctly 
constructed, a surrogate model represents a very close approximation of the behavior 
of the high-fidelity analyses over a given set of parameter ranges. Therefore, design 
space knowledge is significantly improved. 

6.2.2 Chemically Active Cooling 

This study evaluated the use of hydrogen fuel as well as the inert species helium and 
nitrogen for use as coolants to be injected at the wall of the vehicle forebody. Although 
hydrogen would allow a single storage and cryogenic system to be used, its ability to rapidly 
combust presents a practical issue to its implementation. It was shown that at the right 
injection conditions, hydrogen is able to delay ignition significantly, but this benefit may 
have difficulty scaling to very large, orbital-access class flowpaths. One alternative is to use 
helium, which was found to have excellent mitigation properties at the expense of needing 
another storage system. Another expense of helium against hydrogen is that the burning 
of the hydrogen coolant in the combustor may contribute to thrust, whereas helium will 
remain inert throughout the entire flowpath. 

The question then is if there is a coolant that can be used to limit autoignition while 
itself remaining combustible in the combustor? One possibility would be the use of cold 
methane (CH 4 ) at the wall. When methane creates methyl (CH 3 ) species through collision 
with OH and O species, the methyl then has a reaction pathway available to it to recombine 
with monatomic hydrogen to create methane. This reaction takes the form 

CH 3 + H^ CH 4 (88) 

This reaction would be called a “scavenging” reaction — removing hydrogen radicals that 
are important for combustion to proceed and restoring them into methane molecules. By 
introducing methane chemistry into the hydrogen chemistry otherwise occurring at the wall, 
it may be possible to slow down the overall chemistry and delay ignition longer than would 
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be possible with the use of hydrogen alone while still being able to ignite in the combustor. 
Furthermore, less methane might be needed due to its chemical activity versus simply using 
bulk hydrogen to act as a heat sink only. 

Initial pathfinding studies could be limited to the important mechanisms, for example 
the methyl-creating reactions 

CH 4 + OH => CH 3 + H 2 0 (89) 

CH 4 + 0^ CH 3 + OH (90) 

as well as the hydrogen scavenging reaction in Equation 88. It may also be possible to use an 
extinguishing compound such as halon (CFaBr) as the coolant gas, which may only require 
very small amounts to mitigate autoignition compared to hydrogen or methane. The fact 
that it is not combustible, however, and other environmental concerns are detriments to 
using it. 

6.2.3 Effect of Blunt Bodies and Non-Optimal Trajectories 

The two-dimensional notional combustor studied in this work was a purely idealized case 
with no simulation for shock-on-lip conditions or spillage. The actual forebody geometry 
would possibly be subject to shock interactions and a stagnation region near the cowl lip, 
which would be compounded by the fact that the gas may be combustible at that point as 
well. Other possible stagnation regions in the flow would include areas of shock-boundary 
layer interaction as well as engine module walls. These regions may be susceptible to 
instabilities similar to those due premixed, blunt body phenomena observed in numerical 
and experimental literature as well as the present study. 

These phenomena could be studied in the context of the entire flowpath (i.e. a complete, 
integrated flowpath simulation from tip-to-tail) or in more unit problems as was the strategy 
in the current work. For example, the cowl leading edge may be modeled as a wedge of 
variable bluntness, which was studied in a premixed flow to a small degree by Matsuo. 
The effects of both wedge angle and leading edge radius would parameterize the study of 
any instabilities arising in the flow. Module walls — present in the flow whether or not the 
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vehicle is over- or under-sped — separating engine modules could be simulated in a similar 
manner. Boundary layer separation could also be modeled as a blunt body for analysis with 
the quasi-ID framework outlined in this dissertation. 

6.2.4 Integrated Analysis using Variable Geometry 

It is commonly acknowledged that a variable geometry inlet is necessary for a successful 
transatmospheric flight concept at the cost of greater mechanical complexity. The current 
study assumed a fixed inlet geometry designed for a Mach 12 point on a 1500 psf trajectory, 
duly penalizing geometries that spill more air or fuel than others. A variable geometry 
inlet can mitigate some of the issues of fuel-air spillage by actuating the cowl lip to respond 
undersped 2 conditions a flight vehicle experiences as it travels on its design trajectory. By 
moving the cowl closure location, air and fuel that would otherwise be spilled can be directed 
internally for combustion. 

It was observed through the course of this study that, in addition to the injector geome- 
tries themselves, the injected gas does have some contribution to the spillage experienced 
in any particular case. It is hypothesized that varying the primary injection conditions to 
respond to changing fueling requirements as well as varying secondary injection conditions 
to respond to changing cooling requirements throughout a trajectory will have a substantial 
effect on the spillage experienced by a forebody. Therefore, a variable-geometry cowl would 
also need to account for not only the variation in freestream conditions but the conditions 
of the injectants supplied to the forebody. 

A variable geometry nozzle may also be beneficial to the operation of a vehicle. A vari- 
able geometry nozzle, besides altering the shape of the nozzle itself, could vary the location 
of the combustor-nozzle interface in order to begin expanding the flow at different locations 
depending on the trajectory point of the vehicle. This would benefit the thrust generation 
of the vehicle due to the fact that for some arbitrary length combustor, there is a critical 
location where the thrust gains from combustion are outbalanced by the growing viscous 

2 An undersped condition is one where the Mach number is lower than the required Mach number to cause 
a shock-on-lip condition that captures all of the forebody air and ideally causes shock cancellation to occur 
on the forebody shoulder. 
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losses generated at the walls of the combustor. Beginning the nozzle expansion process at 
this “optimum” point would ensure that the greatest amount of thrust is generated. 

6.2.5 Further Development and Application of Quasi-ID Framework 

A significant result of the present research was the development of a tool (QUIVR) for 
the study of blunt body flows of a general purpose nature. Although the current QUIVR 
application is the study of premixed, shock-induced combustion of premixed hydrogen- 
air, the code itself can simulate the stagnation region of any circular blunted-body with or 
without chemistry. Furthermore, the chemistry is of a general nature to allow the simulation 
of different chemical mechanisms. Therefore, the utility of QUIVR extends to other areas 
of research. 

One example research area would be the study of receptivity, which analyses how acous- 
tic disturbances in the main stream interact with boundary layers to cause their transition 
to turbulence. [125]. These acoustic disturbances are commonly found in wind tunnel 
testing, and their understanding is important in order to provide computational results 
that reproduce those in experiment. Furthermore, while the flight environment is typically 
considered to be “quiet” in terms of the generation of acoustic disturbance, facets of the 
flight vehicle (for example, the engine) could source some acoustic disturbance for down- 
stream components. QUIVR could be augmented with low-dissipation schemes to study 
the receptivity phenomenon rapidly. 

6.3 Avenues of Future Experimental Work 

6.3.1 Hypervelocity Validation and Calibration Database 

The experimental simulation of high-enthalpy, hypervelocity flows is of immense perfor- 
mance for the purpose of providing a means of validation and calibration for numerical 
simulations at similar conditions. It was shown here and in previous work that, for exam- 
ple, the selection of turbulent Schmidt number used in RANS simulation is important to 
reproducing the correct physics and performance of geometries and flow phenomena. The 
purpose of this work was not to provide a recommendation on the correct turbulence param- 
eters to use in RANS turbulence models, but auxiliary studies accompanying the body of 
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work presented here show that the simulations are sensitive to the correct selection of these 
parameters. This is generally true, not just for cases of forebody fuel injection, mixing, and 
autoignition as studied here [126]. 

For the purposes of forebody mixing and igniting flows, experimental studies of interest 
would include the testing of generic transverse, ramp, and strut geometries and tracking 
the mixing performance using inert gases or combustible ones, if desired. The simulation of 
film cooling using slots or angled injectors would also yield useful information on the proper 
turbulence modeling of those features in high enthalpy, hypervelocity flows. 

6.3.2 Testing and Evaluation of Forebody Mixing Concepts 

Beyond conducting experiments for the purposes of validating simulations of hypervelocity 
forebody fuel injection, mixing, and autoignition, further experimental attention should be 
given to injectors optimized through a numerical framework like one described in Section 
6.2.1. Experimentally evaluating high-potential geometries yields the maximum benefits 
of both numerical and experimental capabilities by identifying high performance concepts 
through numerical analysis and then verifying that the performance observed actual occurs 
through experimental testing. Such testing would be a necessary step in proof-of-concept 
studies toward the development of a flight test vehicle utilizing forebody injection and 
mixing in the hypervelocity regime of flight. 

Initial experiments on hypervelocity injection and mixing may be conducted on an open- 
plate in a high-enthalpy test facility. Test conditions may reproduce forebody Mach numbers 
while using enthalpies representative of hypervelocity flight. The use of Schlieren imaging or 
a laser-based method such as NO-PLIF may be used in order to observe the flow structures 
in the flow. Initial studies could use a pure helium or helium-air simulant to conduct non- 
reacting mixing studies. A traversing gas sampling system could be used to determine the 
concentration of injectant over a grid of cross-sectional points in order to arrive at estimates 
of mixing efficiency. 
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APPENDIX A 


TAYLOR SERIES EXPANSION OF FLOW VARIABLES 


In order to arrive at a quasi-ID formulation of the Navier Stokes equations, it is necessary to 
recast the flow variables in a form that allows for the expression of a system of equations that 
are a function of distance along the stagnation streamline only. The first step is to create a 
Taylor series expansion of each flow variable about the stagnation (s = 0) streamline. Refer 
to Figure 126 for a schematic of the flowfield. The Taylor series expansion of a generic flow 
variable has the form 


x , ^ x , ^ , d(j)(n, 0) <9 2 </>(n, 0) 2 

0(n, s) = <p(n, 0) H si ^ — s + 


(91) 


ds “ ' d s 2 

where 0 can be any flow variable (e.g. density). Because we are in a nondimensionalized 
system, s = a. The small angle approximation further defines s = sin a. For the sake of 
notational brevity, define 0 l = 0 , 02 = , 03 = '' and 


so on. 


0(n, s) = (n) + 02 (n) sin a + 0 3 (n) sin 2 a + •■■ 


(92) 


A.l Enthalpy, Pressure, and Density 


Because enthalpy, pressure, and density are scalar variables, they are symmetric about the 
stagnation streamline. This has the necessary effect of making odd-ordered terms in the 
expansion equal to zero, thereby reducing the form of these expanded variables to 


H(n,a ) = H\ (n) + H 3 (n) sin 2 a + (n) sin 4 a + ■■ ■ (93) 

p(n,a ) = pi in) + p^ (n) sin 2 a + p$ (n) sin 4 a + • • • (94) 

p(n, a ) = p(n)i + p(n ) 3 sin 2 a + p(n ) 5 sin 1 a + • • • (95) 

A. 2 Velocity 


The components of the velocity normal and tangential to the stagnation streamline behave 
differently in terms of symmetry and therefore needs to be treated separately. The normal 


200 


(u) velocity component is not symmetric in terms of sign about the stagnation streamline. 
As such, even-ordered terms are set to zero in the expansion. The first term (without a 
derivative) is also set to zero because the normal velocity must necessarily be zero along 
the stagnation streamline (ct = 0 ). Thefore, the normal velocity component has the form 

u(n , a) = 112(0) since + u±(o) sin 3 a + uq(o) sin 5 a + ■■ ■ ( 96 ) 

The tangential (v) velocity, unlike the normal velocity, is symmetric about the stagnation 
streamline, thereby eliminating the odd-ordered terms in the expansion. The small angle 
approximation cos(ce) ~ 1 is also used to multiply each term by cos(ce). Therefore, the 
tangential velocity component has the expanded form 

v(o, a) = vi(n) cos a + V3 (n) cos a sin 2 a + v§ (n) cos a sin 4 a + • • • ( 97 ) 

This approximation allows for the successful elimination of a from the reformulated Navier 
Stokes equations and remains exact for a = 0 . 

A. 3 Variable Summary 

The expanded variable set is reproduced here with changed notation such that, for simplicity, 
sequential terms are also sequential in their numbering. 


H(o,a ) = H\(n) + H2 (n) sin 2 a + H3 (n) sin 1 a + • • • ( 98 ) 

p(n, a) = pi(n) + P2(o) sin 2 a + P3(n) sin 4 a + • • • ( 99 ) 

p(n, a) = pi(n) + p2(o) sin 2 a + p3(n) sin 4 a + ■ ■ • (100) 

u(n , a) = u\{n) sinct + 112(0) sin 3 a + 113(0) sin 5 a + • • • (101) 

v(n, a) = v\(n) cos a + V2 (n) cos a sin 2 a + V3 (n) cos a sin 4 a + • • • (102) 


Note that u\ = ^ and P2 = jjji ■ This is the origin of the term “quasi-ID” when referring 
to this formulation, as these terms play an important role in the recast Navier-Stokes 
equations. 
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APPENDIX B 


GRID METRICS AND COMPUTATIONAL DOMAIN 


B.l Computing the Grid Metrics and Grid Jacobian 


Converting the modified Navier-Stokes equations to the computational domain is helpful 
when considering nonuniform grid spacing or for grids that change in time, for example 
in shock fitting schemes. Deriving the grid metrics and grid Jacobian for the quasi-ID 
equations is relatively straightforward as the only variables that must be converted are t (to 
r) and n (to £). The the grid metrics for converting the physical to computational domain 
can be written as 


dr = dt 

d£ = indn + £ t dt 


( 103 ) 

( 104 ) 


where = d^/dn and = d^/dt. Writing these equations in matrix form yields 


1 0 


dt 


dr 



dn 


1 

1 


Similarly, for converting from the computational domain to the physical domain, 


( 105 ) 


where n £ = dn/d £ and n T 


dt = dr 

dn = n^d^ + n T dr 

dn/dr. Writing these equations in matrix form yields 


( 106 ) 

( 107 ) 


1 0 


dr 


dt 

n T 


da 


dn 


This equation can be solved using Cramer’s rule to yield 


( 108 ) 


dr = J 


dt 

dn 


0 


Jn^dt 


( 109 ) 
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di = J 

where the grid Jacobian is given by 

J = 


1 dt 
n T dn 


= Jdn — Jn T dt 


( 110 ) 


1 


1 


(111) 


1 0 

n T 

Equations 109 and 110 can be compared back to Equations 103 and 104 in order to 
arrive at the relations 

Jn^ = 1 

in = J (H2) 

tt — J TL r 

Note that the grid metrics and grid Jacobian are computed at each grid point. 

B.2 Governing Equations in the Computational Domain 

One the grid metrics and grid Jacobian have been computed, it is now possible to convert 
the modified Navier-Stokes equations to the computational domain. We begin with the 
form of the equations 

rt\W r)A/fP 

(113) 


SU» + «M! +qr = 0 


dt dn 

where the superscript p indicates that the vectors are in the physical domain. The partial 
derivatives with respect to t and n can be related to the grid metrics by the relations 

d , ^ — c 9 t d 


dt — & d£ + ^Ot - & ac + dr 

m = Znlk + = 


(114) 


Re-expressing Equation 113 with these relations and substituting Equation 112 yields 


dW 1 U p 3M p Q p 

— ^t - XT 1 7 XT 1 XT I V — 0 


(115) 


dt J dr di J 
Setting U c = U p / J and Q c = Q p / J, where the superscript c indicates a vector in the 
computational domain, yields 


<9U C <9M p dlJP 

+ -^r- - n T ^r + Q 6 = 0 


dr dt 


dt 


(116) 
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It can be shown that, using the product rule, 


51 P dU p n T TT „5r/ 5U p n r 

""'“aT = — ar + u >f = — ar 

which, when setting M c = M — n T U p , finally yields 


5U C 


+ 


5M C 

~w 


+ Q c = 0 


(117) 


(118) 
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APPENDIX C 


IDEAL FUEL-AIR CONDITIONS ON A HYPERVELOCITY 

FOREBODY 


The procedure by which a trajectory is defined is relatively straightforward at the top level 
and is summarized as follows: 

• Define a range of Mach numbers and a trajectory dynamic pressure. 

• Using an atmospheric model at each defined Mach number, find the corresponding 
altitudes whose conditions give the desired trajectory dynamic pressure. 

• Model the effect of a vehicle forebody by processing the freestream flow through two 
subsequent turning angles. 

• Using the conditions on the second vehicle ramp, mix the air stream with fuel in an 
ideal process that conserves mass, momentum, and energy in order to obtain a mixed 
gas state. 

A sample hypervelocity trajectory ranging from Mach 10 to Mach 20 and obtained using 
the previously described process is shown in Figure 141. Hypersonic, airbreathing propulsion 
typically requires trajectory dynamic pressures ranging from 1000 to 2000 psf. A constant 
trajectory dynamic pressure of 1500 psf is chosen for this study, although a varying dynamic 
pressure schedule could alternatively be defined. The low end of the trajectory represents 
the lower limit Mach number where a engine using PMSIC would start operating while 
the upper limit is chosen arbitrarily. The actual Mach number where staging occurs to a 
self-contained rocket to reach orbit is left open. 

Once the the freestream conditions are shocked through a 5° followed by a 6° turning, 
the airflow is mixed with hydrogen fuel while conserving mass, momentum, and energy. 
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Figure 141: A hypervelocity trajectory at a constant dynamic pressure of 1500 psf is shown. 


The three equations to be solved are summarized respectively as 

N 

^ ' PiUiAi = p e u e A e 


N 


22 Ai(Pi + Piuf) = A e (p e + p e u 2 e ) 

i 

N N 

22 PiUiAiCpiTi = (22 PiUiAi)Cp e T e 


(119) 

( 120 ) 
( 121 ) 


where the streamtube cross-sectional areas in Equation 120 are chosen such that the ratio of 
hydrogen mass flow to the mass flow of oxygen in the air is equal to 0.126 for stoichiometric 
flow. The value of the exit area is left to the analyst, but for the purposes of this study it 
was simply set to the cross-sectional area of the air streamtube. The specific heats in the 
energy equation are also assumed to be constant. The final compatibility equation given by 
the ideal gas law 

p = pRT (122) 

The gas properties before and after mixing at the Mach 12 point on a 1500 psf trajectory 
are summarized in Table 24 and are chosen to be practical for an actual airbreathing system. 
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Table 24: An example summary of the pre- and post-mixed properties are shown for 
stoichiometric hydrogen-air flow at the Mach 12 point on a q = 1500 psf trajectory. 



Air 

Hydrogen 

Hydrogen- Air 

A [nr 2 ] 

0.002125 

0.000395 

0.002125 

p [Pa] 

9111 

9111 

9842 

T[ K] 

541 

390 

496 

u [m/s] 

3576 

5815 

3634 


The areas are chosen so that the final mixture is stoichiometric. 
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APPENDIX D 


FINITE-RATE CHEMISTRY MODELS 


Jachimowski (kf = AT n e E / RT ) 



Equation 

A [mol/cm 3 / s] 

Mi] 

E [cal/mol] 

(1) 

H2 + 02 <=> OH + OH 

1.70E13 

0.00 

48000 

(2) 

H + 02 <=> OH + O 

2.60E14 

0.00 

16800 

(3) 

O + H2 <=> OH + H 

1.80E10 

1.00 

8900 

(4) 

OH + H2 <=> H20 + H 

2.20E13 

0.00 

5150 

(5) 

OH + OH <=> H20 + O 

6.30E12 

0.00 

1090 

(6) 

H + OH + M <=> H20 + M 

2.20E22 

-2.00 

0.00 

(7) 

H + H + M <=> H2 + M 

6.40E17 

-1.00 

0.00 

(8) 

H + O + M <=> OH + M 

6.00E16 

-0.60 

0.00 

(9) 

H + 02 + M <=> H02 + M 

2.10E15 

0.00 

-1000 

(10) 

H02 + H <=> H2 + 02 

1.30E13 

0.00 

0.00 

(ID 

H02 + H <=> OH + OH 

1.40E14 

0.00 

1080 

(12) 

H02 + H <=> H20 + 0 

1.00E13 

0.00 

1080 

(13) 

H02 + 0 <=> 02 + OH 

1.50E13 

0.00 

950 

(14) 

H02 + OH <=> H20 + 02 

8.00E12 

0.00 

0.00 

(15) 

H02 + H02 <=> H202 + 02 

2.00E12 

0.00 

0.00 

(16) 

H + H202 <=> H2 + H02 

1.40E12 

0.00 

3600 

(17) 

0 + H202 <=> OH + H02 

1.40E13 

0.00 

6400 

(18) 

OH + H202 <=> H20 + H02 

6.10E12 

0.00 

1430 

(19) 

M + H202 <=> OH + OH + M 

1.20E17 

0.00 

45500 

(20) 

0 + 0 + M <=> 02 + M 

6.00E17 

0.00 

-1800 

(21) 

N + N + M <=> N2 + M 

2.80E17 

-0.75 

0.00 

(22) 

N + 02 <=> NO + 0 

6.40E09 

1.00 

6300 

(23) 

N + NO <=> N2 + 0 

1.60E13 

0.00 

0.00 

(24) 

N + OH <=> NO + H 

6.30E11 

0.50 

0.00 

(25) 

H + NO + M <=> HNO + M 

5.40E15 

0.00 

-600 

(26) 

H + HNO <=> NO + H2 

4.80E12 

0.00 

0.00 

(27) 

0 + HNO <=> NO + OH 

5.00E11 

0.50 

0.00 

(28) 

OH + HNO <=> NO + H20 

3.60E13 

0.00 

0.00 

(29) 

H02 + HNO <=> NO + H202 

2.00E12 

0.00 

0.00 


continued. . . 
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Jachimowski (continued) 



Equation 

A [mol/cm a /s] 

n [1] 

E [cal/mol] 

(30) 

H02 + NO <=> N02 + OH 

3.40E12 

0.00 

-260 

(31) 

H + N02 <=> NO + OH 

3.50E14 

0.00 

1500 

(32) 

O + N02 <=> NO + 02 

1.00E13 

0.00 

600 

(33) 

M + N02 <=> NO + O + M 

1.16E16 

0.00 

66000 


Third body efficiencies: 

(6) : H20 = 6.0 

(7) : H2 = 2.0, H20 = 6.0 

(8) : H20 = 5.0 

(9) : H2 = 2.0, H20 = 16.0 
(19): H20 = 15.0 
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Rogers and Schexnayder ( kf = AT n e E ! RT ) 



Equation 

A [mol/cm 3 / s] 

Mi] 

E [cal/mol] 

(1) 

M + 02 <=> 0 + 0 + M 

0.72E19 

-1.00 

117908 

(2) 

M + H2 <=> H + H + M 

0.55E19 

-1.00 

103298 

(3) 

M + H20 <=> H + OH + M 

0.52E22 

-1.50 

118000 

(4) 

M + H + 02 <=> H02 + M 

0.23E16 

0.00 

-800 

(5) 

M + N02 <=> NO + 0 + M 

0.11E17 

0.00 

64995 

(6) 

M + NO <=> N + O + M 

0.41E19 

-1.00 

149680 

(7) 

M + O + CO <=> C02 + M 

0.30E15 

0.00 

3000 

(8) 

M + H + NO <=> HNO + M 

0.54E16 

0.00 

-596.1 

(9) 

M + H202 <=> OH + OH + M 

0.12E18 

0.00 

45500 

(10) 

M + OH + NO <=> HN02 + M 

0.80E16 

0.00 

-1987 

(11) 

M + OH + N02<=> HN03 + M 

0.13E17 

0.00 

-2200 

(12) 

M + 03 <=> 02 + 0 + M 

0.13E22 

-2.00 

25433 

(13) 

M + HCO <=> CO + H + M 

0.20E13 

0.50 

27400 

(14) 

M + 0 + H <=> OH + M 

0.71E19 

-1.00 

0 

(15) 

H20 + 0 <=> OH + OH 

0.58E14 

0.00 

18000 

(16) 

H2 + OH <=> H20 + H 

0.20E14 

0.00 

5166 

(17) 

02 + H <=> OH + 0 

0.22E15 

0.00 

16800 

(18) 

H2 + 0 <=> OH + H 

0.75E14 

0.00 

11099 

(19) 

H2 + 02 <=> OH + OH 

0.10E14 

0.00 

43000 

(20) 

H + H02 <=> H2 + 02 

0.24E14 

0.00 

695 

(21) 

H2 + 02 <=> H20 + 0 

0.41E14 

0.00 

50470 

(22) 

H + H02 <=> OH + OH 

0.24E15 

0.00 

1887 

(23) 

H20 + 0 <=> H + H02 

0.58E12 

0.50 

57000 

(24) 

0 + H02 <=> OH + 02 

0.50E14 

0.00 

1000 

(25) 

OH + H02 <=> 02 + H20 

0.30E14 

0.00 

0 

(26) 

H2 + H02 <=> H20 + OH 

0.20E14 

0.00 

25000 

(27) 

H02 + H2 <=> H + H202 

0.73E12 

0.00 

18677 

(28) 

H202 + H <=> OH + H20 

0.32E15 

0.00 

8950 

(29) 

H02 + OH <=> 0 + H202 

0.52E11 

0.50 

21062 

(30) 

H02 + H20 <=> OH + H202 

0.28E14 

0.00 

32785 

(31) 

H02 + H02 <=> H202 + 02 

0.20E13 

0.00 

0 

(32) 

0 + 03 <=> 02 + 02 

0.10E14 

0.00 

4790 

(33) 

03 + NO <=> N02 + 02 

0.54E12 

0.00 

2384 


continued. . . 
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Rogers and Schexnayder (continued) 



Equation 

A [mol/cm 3 / s] 

Mi] 

E [cal/mol] 

(34) 

03 + H <=> OH + 02 

0.70E14 

0.00 

1113 

(35) 

03 + OH <=> 02 + H02 

0.90E12 

0.00 

1987 

(36) 

O + N2 <=> NO + N 

0.50E14 

0.00 

75386 

(37) 

H + NO <=> OH + N 

0.17E15 

0.00 

48681 

(38) 

O + NO <=> 02 + N 

0.15E10 

1.00 

38746 

(39) 

N02 + H <=> NO + OH 

0.35E15 

0.00 

1470 

(40) 

N02 + O <=> NO + 02 

0.10E14 

0.00 

600 

(41) 

N02 + H2 <=> HN02 + H 

0.24E14 

0.00 

29000 

(42) 

H02 + NO <=> N02 + OH 

0.30E13 

0.50 

2400 

(43) 

N02 + H20 <=> HN02 + OH 

0.32E13 

0.00 

43714 

(44) 

N02 + OH <=> HN02 + O 

0.21E13 

0.00 

24996 

(45) 

CO + OH <=> C02 + H 

0.70E12 

0.00 

1987 

(46) 

C02 + O <=> 02 + CO 

0.25E12 

0.50 

55040 

(47) 

H20 + CO <=> HCO + OH 

0.65E14 

0.30 

103026 

(48) 

OH + CO <=> HCO + O 

0.58E13 

0.32 

86295 

(49) 

H2 + CO <=> HCO + H 

0.12E14 

0.29 

88242 

(50) 

H02 + CO <=> C02 + OH 

0.15E15 

0.00 

23645 

(51) 

HNO + H <=> H2 + NO 

0.48E13 

0.00 

0 

(52) 

HNO + OH <=> H20 + NO 

0.36E14 

0.00 

0 

(53) 

NO + CO <=> C02 + N 

0.46E09 

0.50 

23983 

(54) 

N02 + CO <=> NO + C02 

0.10E13 

0.00 

27600 

(55) 

NO + H02 <=> HNO + 02 

0.72E11 

0.50 

10928 

(56) 

HNO + O <=> NO + OH 

0.50E12 

0.50 

0 

(57) 

HN03 + O <=> H02 + N02 

0.10E12 

0.00 

0 

(58) 

H02 + N02 <=> HN02 + 02 

0.20E12 

0.00 

0 

(59) 

HCO + 02 <=> CO + H02 

0.10E12 

0.50 

5400 

(60) 

03 + H02 <=> 02 + 02 + OH 

0.10E12 

0.00 

2800 


Third body efficiencies: 

(1) : 02 = 4.0, 0 = 10.0, H20 = 2.0 

(2) : H2 = 2.0, H = 2.0, H20 = 2.0 

(3) : H20 = 6.0 

(4) : H2 = 2.0, H20 = 13.0 

(8) : H20 = 3.0 

(9) : H20 = 6.0 
(12): 02 = 1.5 
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O Conaire (Using Troe Falloff Function) 



Equation 

A [mol/cm a /s] 

n [1] 

E [cal/mol] 

(1) 

H + 02 <=> O + OH 

1.91E14 

0.00 

16440 

(2) 

O + H2 <=> H + OH 

5.08E04 

2.67 

6292 

(3) 

OH + H2 <=> H + H20 

2.16E08 

1.51 

3430 

(4) 

O + H20 <=> OH + OH 

2.97E06 

2.02 

1340 

(5) 

H2 + M <=> H + H + M 

4.57E19 

1.40 

105100 

(6) 

O + O + M <=> 02 + M 

6.17E15 

0.50 

0.00 

(7) 

O + H + M <=> OH + M 

4.72E18 

1.00 

0.00 

(8) 

H + OH + M <=> H20 + M 

4.50E22 

2.00 

0.00 

(9) 

H + 02 + M <=> H02 + M 

k ffi = 3.48E16 

0.41 

1120 


H + 02 <=> H02 

k f = 1.48E12 

0.60 

0.00 

(10) 

H02 + H <=> H2 + 02 

1.66E13 

0.00 

820 

(11) 

H02 + H <=> OH + OH 

7.08E13 

0.00 

300 

(12) 

H02 + O <=> OH + 02 

3.25E13 

0.00 

0.00 

(13) 

H02 + OH <=> H20 + 02 

2.89E13 

0.00 

-500 

(14a) 

H02 + H02 <=> H202 + 02 

4.20E14 

0.00 

11980 

(14b) 

H02 + H02 <=> H202 + 02 

1.30E11 

0.00 

1629 

(15) 

H202 + M <=> OH + OH + M 

kffi = 1.27E17 

0.00 

45500 


H202 <=> OH + OH 

k f = 2.95E14 

0.00 

48400 

(16) 

H202 + H <=> H20 + OH 

2.41E13 

0.00 

3970 

(17) 

H202 + H <=> H2 + H02 

6.03E13 

0.00 

7950 

(18) 

H202 + O <=> OH + H02 

9.55E06 

2.00 

3970 

(19) 

H202 + OH <=> H20 + H02 

1.00E12 

0.00 

0.00 


H202 + OH <=> H20 + H02 

5.80E14 

0.00 

9560 


Third body efficiencies: 

(5): H20 = 12.0, H2 = 2.5 


(6) : H20 = 12.0, H2 = 2.5, Ar = 0.83, He = 0.83 

(7) : H20 = 12.0, H2 = 2.5, Ar = 0.75, He = 0.75 

(8) : H20 = 12.0, H2 = 0.73, Ar = 0.38, He = 0.38 

(9) : H20 = 14.0, H2 = 1.3, Ar = 0.67, He = 0.67 
(15): H20 = 12.0, H2 = 2.5, Ar = 0.45, He = 0.45 
Troe parameters: 

(9): A = 0.5, T3 = 1.0E-30, T1 = 1.0E30, T2 = 1.0E100 
(15): A = 0.5, T3 = 1.0E-30, T1 = 1.0E30 
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San Diego Hydrogen Submechanism (Using Troe Falloff Function) 



Equation 

A [mol/cm 3 / s] 

Mi] 

E [cal/mol] 

(1) 

H + 02 <=> OH + 0 

3.52E16 

-0.70 

17069.79 

(2) 

H2 + 0 <=> OH + H 

5.06E04 

2.67 

6290.63 

(3) 

H2 + OH <=> H20 + H 

1.17E09 

1.30 

3635.28 

(4) 

H20 + 0 <=> OH + OH 

7.60E00 

3.84 

12779.64 

(5) 

H + H + M <=> H2 + M 

1.30E18 

-1.00 

0.00 

(6) 

H + OH + M <=> H20 + M 

4.00E22 

-2.00 

0.00 

(7) 

O + O + M <=> 02 + M 

6.17E15 

-0.50 

0.00 

(8) 

H + O + M <=> OH + M 

4.71E18 

-1.00 

0.00 

(9) 

O + OH + M <=> H02 + M 

8.00E15 

0.00 

0.00 

(10) 

H + 02 + M <=> H02 + M 

k ffi = 5.75E19 

-1.40 

0.00 


H + 02 <=> H02 

k f = 4.65E12 

0.440 

0.00 

(11) 

H02 + H <=> OH + OH 

7.08E13 

0.00 

294.93 

(12) 

H02 + H <=> H2 + 02 

1.66E13 

0.00 

822.90 

(13) 

H02 + H <=> H20 + O 

3.10E13 

0.00 

1720.84 

(14) 

H02 + O <=> OH + 02 

2.00E13 

0.00 

0.00 

(15) 

H02 + OH <=> H20 + 02 

2.89E13 

0.00 

-497.13 

(16) 

OH + OH + M <=> H202 + M 

fc/,o = 2.30E18 

-0.90 

-1701.72 


OH + OH <=> H202 

k f = 7.40E13 

-0.37 

0.00 

(17) 

H02 + H02 <=> H202 + 02 

3.02E12 

0.00 

1386.23 

(18) 

H202 + H <=> H02 + H2 

2.30E13 

0.00 

7950.05 

(19) 

H202 + H <=> H20 + OH 

1.00E13 

0.00 

3585.09 

(20) 

H202 + OH <=> H20 + H02 

7.08E12 

0.00 

1434.03 

(21) 

H202 + O <=> H02 + OH 

9.63E06 

2.00 

3991.40 


Third body efficiencies: 

(5, 6, 7, 8, 9): H2 = 2.5, H20 = 12.0 
(10): H2 = 2.5, H20 = 16.0 
(16): H2 = 2.0, H20 = 6.0 
Troe parameters: 

(10): A = 0.5, T3 = 1.0E-30, T1 = 1.0E30 
(16): A = 0.735, T3 = 94, T1 = 1756, T2 = 5182 
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